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' Don State Technical University (Rostov-on-Don, Russian Federation) 
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Introduction. When developing ocular prostheses, a number of problems arise, one of which is the construction of the 
connection between the hard optical part and the soft corneal tissue. Their Young's modules can differ by three orders of 
magnitude. In this case, the problem arises of creating an intermediate layer, possibly with gradient properties, whose 
purpose is to exclude injury to soft biological tissues. Two types of keratoprostheses are considered: the first type with a 
support plate and the second type with an intermediate functionally gradient layer. The stress-strain state of the prosthe- 
sis 1s calculated for the first type. For the second type, analytical and finite element modeling of the interaction of a cy- 
lindrical optical prosthesis, an intermediate inhomogeneous layer, and the cornea was carried out in the elastic media. 
Two versions are considered: discounting the curvature (circular plate or plate) and with account of the curvature 
(spherical dome or shell). The work objective is to study the stress-strain state of the keraprosthesis and cornea in the 
contact area. 

Materials and Methods. Mathematical models of the structures under consideration are the boundary value problems of 
the linear elasticity theory. The analytical solution is constructed for a simplified model in the form of a composite cir- 
cular plate. Spatial three-dimensional problems and axisymmetric problems are solved by the finite element method. 
Finite element modeling of the considered structures was performed in the CAE package ANSYS and ACELAN. 
Results. CAD models of keratoprostheses with conditions of fixing and loading are constructed. The load acting on the 
keraprosthesis under the effect of intraocular pressure was determined. The stress-strain state of the keratoprosthesis 
and cornea elements was calculated. Special attention was paid to the area of its contact with the keratoprosthesis. 
Discussion and Conclusions. The results of calculating the axial displacements and mechanical stresses in the first type 
of keratoprosthesis show that the selected geometric parameters meet the kinematic and strength requirements. The pro- 
posed models of the deformed state of soft biological tissues provide assessing their injury when using a keratoprosthe- 
sis of the second type, as well as selecting the geometric parameters and gradient properties of the intermediate layer. 


Keywords: ocular prosthesis, inhomogeneous elastic properties, plate, shell, finite element method. 
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Introduction. A keratoprosthesis (KPro) is a cell-free artificial implant designed so that a cylindrical frame 





holds the optics. Keratoprosthesis replaces the removed part of the cornea. Previously, corneal transplants had high rates 
of infection and rejection. In the late 1980s, the “core-skirt” design (a bio-integrated “skirt” surrounds the optics) was 
the most widespread one. Not only the sizes, but also the location of the pores in the porous skirt were important. Kera- 
toprostheses, such as AlphaCor (formerly known as Chirila), were polymethylmethacrylate devices with a central opti- 
cal region fused to the surrounding porous skirt. 

Modern keratoprostheses consist of an optical element and a supporting plate. The optical transparent element 
has the shape of a cylinder or lens. The supporting plate connected to the optical element can have various shapes: a 
ring with holes, a spoked wheel, ears, or amoeboid legs. The book by S. N. Fedorov [1] describes various types and 
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forms of keratoprostheses and fasteners. The main complication after keratoprosthetics is aseptic necrosis of the cornea 
which develops in front of the supporting attachment of the implant. This complication often causes rejection of the 
keratoprosthesis [1]. The major cause of aseptic necrosis is blocking the flow of vital substances into the layers of the 
cornea, which are located above the keratoprosthesis support. In this regard, for the manufacture of the supporting plate, 
it is required to use a biocompatible material that will allow corneal tissues to grow through the supporting mount. 

American researchers have described a microporous support plate' for keratoprosthesis made of stretched poly- 
tetrafluoroethylene (PTFE). Its structure has a configuration in the form of polymer nodes, which are connected by fi- 
brils 7-8 um long. Too small pores prevent the corneal tissue from ingrowth. Therefore, it is required to additionally 
penetrate the prosthesis with pores perpendicular to the surfaces of the supporting plate with a diameter of 20-150 um 
(preferably 50 wm). The thickness of the supporting plate should be about 0.2 um, but not more than 0.3 wm. 

In the model with a ring-shaped’ supporting attachment, the optical element is made of a transparent substance 
such as polymethylmethacrylate (PMMA), the attachment is made of a hydrophilic porous material with fibrous struc- 
ture through which corneal tissue can ingrow. A supporting plate with a thickness of 0.15—0.30 um is formed from this 
fiber. A ring with an external diameter of about 9.5 mm is placed on the optical element. Its front part is cylindrical, the 
back has a conoid, whose larger section is directed inward. The supporting plate has elastic characteristics required to 
prevent rejection of the keratoprosthesis and necrosis due to the pressure of the surrounding eye tissues. 

In another design’, the supporting plate consists of two parallel shells spaced 0.40.7 um apart. The optimal 
distance is 0.3+0.02 mm. Conical shells are attached to a cylindrical optical element. The material of the supporting 
plates should be bio-populated; therefore, it is assumed that its porosity is 50% or higher, the pores should be open, and 
their diameter is about 20—100 wm. 

Elements of the keratoprosthesis* proposed by Russian scientists (Fig. 1) are as follows: 

— supporting plate — 1 (radius of curvature is 7-10 mm, external diameter is 7—12 mm); 

— optical part — 2 (optical element in the form of a cylinder with rounded ends that serve to select diopters); 

— cylindrical surface of the optical threaded element — 3; 

— connecting washer — 4 (rigidly attached to the supporting plate). 





Fig. 1. Scheme of keratoprosthesis 


The keratoprosthesis is implanted in the intralamellar corneal pocket. The supporting plate is located inside the 
cornea, and the optical element penetrates the entire cornea. 

The implant for strengthening the cornea” is a plate — round, oval, in the form of a shamrock, chamomile, or 
convex-concave lens with a curve radius of 7-10 mm. There may be a cylindrical hole in the center of the lens. In the 
structure of the implant material, the volume void rating is 15—40 %, their specific surface area is 0.25—0.55 mm’/ um’, 
the average distance between voids is 25—50 um, and the mean volume chord is 8—25 um. The implant is a convex- 
concave lens whose curve radius corresponds to the curve radius of the patient's cornea. In the central part, the implant 
thickness is 0.3—0.7 mm, and it is 0.01 mm at the edges. 

In another modification of the keratoprosthesis°, in contrast to the one described above, the supporting block is 
made of porous titanium nickelide, the curvature of the plate coincides with the curvature of the cornea. The supporting 
plate is installed on the surface of the cornea and is retentively fixed with a sclera allograft. The cylindrical optical ele- 
ment is connected to the supporting plate by a tight fit method. This approach enables to treat leukomata more effective- 


' US Patent no. 5713956, M. cl. 7 A 61 F 2/14, 1998. 

* US Patent no. 5489301, M. cl. 6 A 61 F 2/14, 1996. 

> US Patent no. 6106552, M. cl. 7 A 61 F 2/14, 2000. 

* Keratoprosthesis. RF Patent no. 2270643, 2006. (In Russ.) 

> Implant to strengthen the cornea. RF Patent no. 2270642, 2006. (In Russ.) 

° Keratoprosthesis and method of surgical treatment with its help. RF Patent no. 2367379, 2009. (In Russ.) 
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ly, thanks to the rapid and successful integration of the porous implant (base) with the surrounding tissues and mechani- 
cally strong fixation of the optical element, which reduces the risk of keratoprosthesis reposition. 

The developer of the combined keratoprosthesis' notes that the supporting part of the implant should maximal- 
ly resist pushing pressure of the intraocular fluid on the optical part of the keratoprosthesis and minimally deform the 
appropriate tissues. The material of the supporting part should not only be biocompatible, but also structured so as to 
provide the fusion of the above and below layers of the cornea, which are separated during the implantation of the kera- 
toprosthesis. Keratoprostheses with an optical element made of a transparent substance (for example, PMMA) and a 
supporting plate in the form of a ring made of hydrophilic porous material PTFE or polyethylene” were fore-mentioned. 
Their use provides fusion of the above and underneath layers of the cornea due to the high porosity of the material 
(from 50 %) and the pore diameter up to 100 um. The main disadvantages of such structures are excessive flexibility 
and low mechanical strength. These characteristics do not provide reliable retention and fixation of the installed kerato- 
prosthesis with large optics for a long time. 

The optical part of the combined keratoprosthesis is made of a transparent elastic polymer saturated with UV 
adsorbent. In shape, it is a removable bolt with a diameter of 5-6 mm. The front end is fungiform, while the others 
are spherical (aspherical). The haptic part is a supporting disk with an outer diameter of 8-12 mm and a thickness of 
0.3—0.9 mm. Its side surface is threaded. For this part, a porous (perforated) solid (elastic) areactive polymer (metal) 
with a volume void rating of no more than 50 % is used. The front part of the haptic is made of biological materials, the 
back part is made of artificial biocompatible materials. The supporting disk of the back part is flat and connects perpen- 
dicularly to the supporting hollow cylinder located in the center. The connection is provided by a groove and / or flange 
on the outer rear of the cylinder. The profile and diameter of the internal thread of the cylinder correspond to the thread 
of the optical bolt. The support disk can pass into a sleeve that wraps around the supporting cylinder from the outside. 
The supporting plate is connected at the front with a round biodisc, which is cut out of freshly removed or preserved 
sterile biological materials. They are multi-layered combinations of biological tissues. Before joining, parts of the kera- 
toprosthesis are sterilized. The connection is carried out under sterile conditions. Assembly is performed immediately 
before prosthetics or in advance (in this case, the implant 1s stored until the surgery). The external diameter of the bio- 
disc should be greater than the external diameter of the supporting disk by 0.3—1 mm. In the center of the biodisc, a co- 
axial bore is made with a diameter corresponding to the outer diameter of the cover sleeve or the supporting cylinder. 

For penetrating keratoplasty, models are developed using various materials, structures, and surgical methods 
[2]. However, the problem is not completely resolved. This is confirmed by cases of corneal blindness in patients with 
repeated graft failure or severe damage to the eye surface [3—6]. The issues of keratoprosthetics in diseases of the cornea 
and ocular surface are considered in [7]. Review on the application of Boston KPro is presented in [8]. In [9], complica- 
tions with wide consequences in the use of keratoprostheses are discussed. 

In this paper, we consider two types of keratoprostheses with a cylindrical optical element. The diagram of the 
first type is shown in Fig. 1. In the second type, there is some intermediate layer with non-uniform mechanical proper- 
ties on the cylindrical surface. Fig. 2 shows a diagram of the half-axial section of the second-type KPro disregarding the 
curvature (a) and with account for the cornea curve (b). The layer with functionally gradient properties is greyed out. It 
is needed for non-traumatic contact of the optics with the soft tissue of the cornea. 





a) b) 


Fig. 2. KPro half axial section: plate model (a), dome model (b) 


Analytical and numerical models of keratoprosthesis implantation in the cornea are constructed using the finite 
element method. The stress-strain state of the cornea in the contact area is studied. This defines the objectives of this 
work. For the first-type KPro, they include investigations on its stress-strain state under the action of intraocular pres- 
sure, assuming that the outer edge of the supporting plate is fixed. The objectives for the second-type KPro are to study 


' Combined keratoprosthesis. RF Patent no. 2707646, 2019. (In Russ.) 
* US Patents no.5489301, no. 6106552. 
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the stress-strain state of the cornea in the vicinity of the contact with the keratoprosthesis, whose outer layer is a struc- 
ture that reduces soft tissue injury. This intermediate (interface) layer is modeled as a hollow cylinder with functionally 
gradient mechanical properties. 

Three models are considered in the paper. 

I. Analytical — for an implanted prosthesis based on the bending of a composite circular plate. 

II. Finite element — for an implanted prosthesis based on a composite circular plate. 

III. Finite element model — for an implanted prosthesis based on a composite spherical dome. 

Modeling problems are solved for: 

— keratoprosthesis of the first type (A), 

— keratoprosthesis of the second type (B). 

Materials and Methods. The continual formulation of the tasks. The general mathematical formulation of the 
problem under study (for problems A and B) is a static boundary value problem of elasticity theory for a composite iso- 
tropic body. 

For the unknown u = Cee U;) of the displacement vector component, the system of differential equations 


has the form [10]: 


1 
G, , =9, && = rae +u,,); (1) 


where o.. 


i> Ci» &y ae components of stress tensors, elastic constants, and deformations, respectively. 


Boundary conditions are set for the displacement and stress vectors on the corresponding surfaces S, and S; : 
u. | =u; (X02), xesS, (2) 
(; Is, = OM, Is =9(X,t), xeES,, (3) 
where n, are coordinates of the outer normal unit vector. 


In addition, for problem B (see Fig. 2) subbodies have different properties, namely: 


— two of them (cylindrical optical prosthesis and cornea) are homogeneous, with elastic moduli Ev, and E7,; 


— functionally gradient elasticity modulus of the interface layer — Er, = Er,(r) . 


The right side is fixed, and the left side has symmetry conditions. The lower boundary is affected by uniform 
pressure, which corresponds to excessive intraocular pressure compared to atmospheric pressure. On the interface 
boundaries, continuity conditions are set. 

To study the model (I), the boundary value problem (1) — (3) is reduced to a system of ordinary differential 
equations (4), (5), and (6) for the first, second, and third sections, respectively, relative to the angle of rotation of the 
normal with a power-law dependence on the radius of the elastic modulus of the second section [11]: 


at 24 (4 _ _i 3 
(s a(7)) ré+ (< a(r)] r—O(r) = > Klar ; (4) 
d? d d 1 K2q(r?-rZ)r?2—4% 
(s a(7)) r-+a (< a(r)] r+ad,0(r) + (= a(r) | r—O(r) = 1 Real (5) 
2 0. aee (= a(r)] p06) = 6G? ar (6) 
dr2 dr 2 a eee? 
where K =1/D, D is cylindrical stiffness. 
Research Results 
Analytical solution 
General solution to the system (4) — (6) is presented below. 
The first section: 
O(r) =rC, + —Klqr? + = 
The third section: 
6(r) = < +rc,—- — K3q(4 In(r)r? — r?)r. (7) 


The second section: 
Q(r) = pozatgV a? 4082 +4 Cy + pozatz la? = 4082 +4 C, + 
1 K,qr+~% (((r? —77)9, —r? + 3rl”)a— sri?) 


2 (8 + (9, — 3)a)a(9, — 1) 
at 


Mechanics 


353 


Advanced Engineering Research 2020. V. 20, no. 4. P. 350-359. ISSN 2687-1653 





Ery 
Er2 h3 _ In (ry) 


_ oa 
oe ae @ a 12(1-02)? * ~ In) 


Arbitrary constants C, of the general solution (7) are determined from the boundary conditions at the right end 





and the blending conditions, and C,, = 0 for the first section. 


Finite element solution 

Numerical simulation was performed in finite element packages ANSYS and ACELAN [12, 13]. 

Numerical results for KPro of the first type. Fig. 3 shows a CAD model of the keratoprosthesis with the 
conditions of fixing and loading. 





‘s 
0,000 9,000 (rim) 
———— | 


4/500 


Fig. 3. Geometric model with fixing and loading conditions 


The outer edge of the supporting plate is fixed, and the optical element is affected by distributed pressure. 
Normal intraocular pressure is between 12 and 21 mm/Hg. The pressure of 45 mm/Hg 1s critical. In this case, an opera- 
tion can still save the eye. With higher pressure, it is almost impossible to save the eye. The current keratoprosthesis has 
actual load of 46 mm/Hg. 

Fig. 4 shows the distribution of the axial displacement (a) and von Mises stress (b) on the finite element grid of 
the model. 





0,000 10,000 (rm) Ls : 
a 





a) b) 


Fig. 4. Distribution of the stress-strain state characteristics on the finite element model: axial displacement (a) 
and stresses according to von Mises (5) 
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Maximum stresses occur at the junction of the optical element and the supporting plate, but their values do not 
exceed the strength limit of the selected material. 

Next, we consider the calculation results for problem B in the axisymmetric formulation. So, Fig. 5—7 show the 
components of the stress-strain state (II). 
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a) b) 


Fig. 5. Calculation of displacements: finite element grid on a deformed structure (a), 
distribution of the displacement vector (5) 
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Fig. 6. Calculation of stress tensor components. Dependence on radius: radial stresses at the upper boundary (a), at the interface 


boundary with the cornea (b), angular stresses at the upper boundary (c), at the interface boundary with the cornea (d) N/mm* 
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Fig. 7. Distribution of components of stresses inside the area: radial stresses (a), 
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The stress-strain state characteristics of model III are shown in Fig. 8—11 
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Fig. 9. Calculation of radial stresses: distribution in the region (a), at the interface 
border with the cornea (b) 
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Fig. 10. Calculation of angular stresses: distribution in the region (a), at the interface border with the cornea (5) 
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Fig. 11. Calculation of shear stresses: distribution in the region (a), at the interface border with the cornea (5) 


Discussion and Conclusions. Two types of keratoprostheses are considered, and their mechanical and mathe- 
matical models are constructed. For the first type of KPro with a supporting plate, its stress-strain state was studied, and 
a finite element model was created in ANSYS. It is shown that at maximum eye pressure, the resulting maximum 
stresses do not exceed the strength limits of the selected materials. For the second type of KPro, an analytical solution 
disregarding curvature and a finite element solution in the ACELAN package with corneal curve are constructed. The 
characteristics of the stress-strain state, including those at the interface border with the cornea, are calculated. This pro- 
vides evaluating its injury and selecting geometric parameters and gradient properties of the intermediate layer. 
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Frequency function and damping function in assessment of dynamic satis 
processes in mechanical oscillatory systems with symmetry 





A. V. Eliseev Oe 
Irkutsk State Railway Transport University (Irkutsk, Russian Federation) 


Introduction. A new approach to the formation of the methodological basis of system analysis in the application to the 
problems on mechanical oscillatory structure dynamics is considered. The study objective is to develop a method for 
evaluating properties of the mechanical oscillatory systems with account for viscous friction forces based on frequency 
functions and damping functions that depend on the so-called coefficient of connection forms, which is the ratio of 
characteristics of generalized coordinates. 

Materials and Methods. The graphoanalytical methods used for evaluating the dynamic properties of mechanical 
oscillatory two-degree-of-freedom systems are based on determining the extreme values of the frequency functions and 
the damping function, which are determined from the relations between the kinetic, potential energy and the values of 
the energy dissipation function. Mathematical models are based on Lagrange formalism, matrix methods, and elements 
of the theory of functions of a complex variable. 

Results. A method is proposed for constructing frequency functions and damping functions for a class of mechanical 
oscillatory two-degree-of-freedom systems based on the analytical expressions that reflect features of the ratio of the 
potential and kinetic energy of the system considering viscous friction forces represented by the dissipative function. 
General analytical expressions for the frequency function and the damping function are derived. Graphoanalytical 
analysis of extreme properties of the corresponding frequency functions and damping functions is performed for 
mechanical vibrational systems with elastic-damping elements with symmetry properties. The results of numerical 
experiments are presented. A criterion for classifying frequency functions and damping functions based on the 
topological features of the graphs of the corresponding functions is proposed. 

Discussion and Conclusions. The developed method for constructing frequency functions and damping functions can be 
used to display the dynamic features of mechanical oscillatory systems. The proposed matrix method for constructing a 
frequency-damping function for a two-degree-of-freedom system can be extended to the mechanical vibrational systems 
considered in different coordinate systems. 


Keywords: mechanical system, dynamic connections, frequency function, damping function, connectivity of movement, 
extreme properties, oscillation, viscous friction. 
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Introduction. Considerable attention is paid to the methods of using mechanical vibration systems as 
calculation schemes in the problems on evaluating the dynamic properties of technical objects operating under intense 
vibration loads [1-9]. Methods based on extreme properties of the ratio of potential and kinetic energy can be 
considered as common approaches to evaluating the dynamic properties of mechanical oscillatory systems [10, 11]. 
Methods based on energy relations have been developed in the use of the frequency function as a function of the 
coefficient of connection forms of the mechanical system coordinates to evaluate the dynamics of mechanical 
oscillatory systems disregarding friction forces [12—15]. 

At the same time, methods for evaluating the dynamic properties of mechanical oscillatory systems, with 
account for the viscous friction forces based on the frequency function, require detailed representations depending on 
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the viscous friction value. In particular, this is due to the fact that for systems with aperiodic motion, the concept of 
oscillation frequency may lose its meaning. 

The proposed work is devoted to the creation of a method for evaluating the properties of mechanical 
movements based on the development of the concept of the frequency function, when an additional damping function 
reflecting the features considering viscous friction forces depending on the coefficients of connection forms, is 
introduced. 


Materials and Methods. Free motions of a mechanical elastic-dissipative system with concentrated two- 
degree-of-freedom parameters are considered. The schematic diagram of the system is shown in Fig. 1. 


v7] Y2 


b, bo b, 
My] MN 


fy ome ome ky 


Fig. 1. Mechanical oscillating system with account for viscous friction 


Generalized coordinates y,, y, denote the displacement of mass-inertial elements m,,m, relative to the static 


equilibrium position. Kinetic energy 7 , potential energy II and scattering function F’ have the form: 


1, 1, 
T= Sm +m yy (1) 
1 2 1 2 1 2 
=e 5 22 eho a Mi) : (2) 
F -_ 1 - 2 4 1 b . . 2 l, - 2 3 
~ > V1 3 0(¥2 - Vi) 5 a2. + (3) 
The system of Lagrange equations of the second kind has the form: 
d oT Ol OF 
—— +— +— =0); 
Ue, OG is 
d oT OW OF _ 


——— +——_ + = Q,, 
dt Oy, Oy, Op, 
After substituting expressions T ,I1 , F', the system of differential equations (4) takes the form: 
ms + (By +B), — Oy + Ky +): — Ko. = 9; (5) 
My Vy + (By + By )V_ — Oy Y, + (Ky + ky). — ky, = 9. 
The forms of free motions of the presented system (5) are generally determined by the eigenvalues of the 
characteristic equation, taking into account their multiplicity. The case of simple roots is considered. Thus, let the 
solution y, = y,(t), ¥, = y(t) of the system (5) be represented as: 


— 
— 


y=Ye", (6) 


~ |M|. , Meee 
where y= is the solution vector, Y = 
y4 2 


is the numeric vector, p=o+ j@ is the complex parameter, ¢ is the 
time variable. It is assumed that the initial conditions are consistent with the type of solution you are looking for (6). 

The task is to construct and evaluate the extreme properties of functions that display the characteristics of the 
proper motions of the system with account for the viscous friction forces. 

Research Results 

1. Construction of the frequency function and the dissipation function based on the energy ratio. The 
system (5) in the notation (6) has the form: 


mp +(by +b) ptk, +k, —by p — ky 
—b, p-k, m,p +(b, +b,)ptk, +k, 
We introduce the notations: 


ae 7 
[=o (7) 
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|" ae —b, Je=feg —k, | (8) 
0 m, —bh, 6, +6, —k,  k, +k, 
Considering (8), the matrix relation (5) takes the form: 

(p’A+ pB+C)Y =0. (9) 


On the basis of the matrix relation (9), various scalar equations can be obtained. With account for their extreme 
properties, in turn, the properties of the solution p of the equation (9) can be determined. In particular, scalar 


multiplication of the left and right parts of the equation (9) by vector Y results in the scalar expression: 
py (AY,¥)+ p(BY,Y)+(CY,Y}=0. (10) 
Using substitution _p =o + jw, we can write (10) as: 
(co? -@ + 2 jow)( AY, Y)+ (o+ jo) (BY, Y)+ (ee Y) =0, (11) 
Let the following relation be fulfilled for the coordinates of vector Y : 
Y, =a, (12) 
where o is the coefficient of the connection form. In this case, vector Y can be presented as: 


Vera. (13) 
l 

where @ -| is vector defined by the connectivity coefficient a. After substituting (13), the expression (11) takes 
Oo 


the form: 
(o -@ +2jow)A, +(o+ jo)B,+C, =0, (14) 
where A, = (AG, &) , B= (Ba, a) B= (Ba, a) are scalar functions of the coefficient of the connection form a. The 
equation (14) can be presented in the equivalent form: 
0 A, =0°A,+0B,+C,; 
oe +oB, =0. 


We find the solution to the system (15) with respect to o,@, as functions a. Function @(a) is considered as a 


(15) 


frequency function, o(a) 1s the damping function. Features of the solution to the system are determined by the sign of 
the expression o°4, +oB,+C,. 

Considering dissipation. The level of energy dissipation can be characterized by discriminant B,* -44_C,. 
Under the following condition: 


B’<4A,C (16) 


a~a? 


which is understood as a small amount of friction, the solution (15) can be presented as: 


2 Cc. B, 
o =—~-|—* | ; 
A. \2A. 
B, 
2A, 
It should be noted that when the conditions of smallness of the viscous friction forces (16) are met, the equation is 
performed: 


(17) 

















Cc 
2 2 
wo +o0°=—. 18 
A. (18) 
Under the condition of “large viscous friction forces”: 
B,* >44,C,, (19) 
the solution (15) can be presented as: 
2 2 
B B C B B C 
o=0, (20) o,(a) = -—~- — | ——* ,o,(a) =-—* + —| ——“, 21 
29) oe) 2A, Ga A, 20) 24, Ga A, e) 
The damping function has two components o,(a@), o,(a) , such that: 
B C 
Oa cr alae (22) 
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Thus, depending on the level of viscous friction forces, different types of presentation of the frequency 


function and the damping function are possible. If B.*-4A4,C, <0, then the frequency function @°(a) and one 
component of the damping function o(a) are defined. If B,*-44,C, >0, then it is assumed that the frequency 
function @*(a) takes zero values, and the damping function has two different negative components o,(a), 6,(a@). 


As for the condition B,* —44,C, =0, it requires a separate analysis. However, the condition B,* -44,C, =0 


can be interpreted as a boundary between two different modes of motion of a mechanical system. 

The presented analytical expressions of the frequency function, the damping function, and the conditions of 
“small” and “large” viscous friction forces can be detailed when considering specific options of mechanical oscillatory 
systems obtained on the basis of a two-degree-of-freedom system. 

2. Frequency function and damping function for a mechanical two-degree-of-freedom system. Parameter 
options for the mechanical system shown in Fig.1 are considered. It is supposed that a set of boundary parameters 
separating the modes of motion for small and large forces of viscous friction is determined from the equation: 


B =44.C,, (23) 
where: 
a ue. (24) 
B =(b +b, )a° -2ab, +b, +5,, (25) 
C .=(k, +k, )a° -2ak, +k, +k,. (26) 
The conditions of smallness of viscous friction forces have the form: 
B,” <44,C,. (27) 


In this case, frequency function (a) and damping function o(a): 


2 
o°(a) = (k, +k, )o? Deas +h +k, (34 +b, )o’ ~2aby +b, “4 
m, +m,o 2 m, +m, (28) 
2 
a 1G + b, )o es +5, 
2 m, +m,O 


The conditions of large viscous friction forces have the form: 
B,* >44,C,. (29) 


Under the conditions (29), functions @* and o(a) have the form: 


1G +5, )a* —20b +h +b, _ 


Z m, +m,a° 


0, (a) =F 





10 +b)0? -20b, +45) Ky tho? 20k, thy th, a 
2 m, +m, m, +m,a” 
1 (b, +b, )a* —2ab, +b, +b 
5, (a) =—-— 0 2 . 0 ae 
Z m, +m,o 





2 
1 (by +b, a? —2ab, +b, a (k, +k,)a°? —2ak, +k, +k, 
a tage «|  . “mec iao. ° 
Z m, +m, m, +m,o 
The presented expressions reflect motions in the form of exponential decrease in the absence of fluctuations. 
3. Features of frequency functions and damping functions for symmetric mechanical oscillatory systems. 
We consider a mechanical oscillating system with elastic damping elements, whose parameter values are imposed by 


symmetry conditions in the form 6, =b, =b, =b, k, =k, =k, =k . The schematic diagram is shown in Fig. 2. 
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Fig. 2. “Symmetrical” mechanical system 


The system of differential equations (5) has the form: 


es + 2b}, — by, + 2ky, —ky, = 0: - 
mii, + 2b, —bj, + 2ky, —ky, = 0. 
Functions A,, B,, C,, can be presented by the expressions: 
A .=m+mo’, (32) 
B .=2b(a’ -a+l), (33) 
C .=2k(a°-a+l). (34) 


On the basis of the presented components, frequency function and damping function can be constructed, and 
the conditions for the smallness of the viscous friction forces can be formulated. 

Accounting for viscous friction forces. The condition of smallness of the friction forces can be presented from 
the inequality: 


2 

B 

—|< Si , (35) 

2A, A, 
After substituting the functions (32) — (34), the condition of smallness of the friction forces (35) can be written 

as: 
Yo <M(a), (36) 
: 1 om, +m,a° 


where y, = ae is a generalized viscoelastic parameter, M(a) = is a generalized mass-inertia coefficient 


2 a -a+l 
that depends on the shape coefficient a. The graph of function M, for each fixed y, defines a set of values for a, at 


which the condition of smallness of the friction forces is satisfied. 
As an example, Fig. 3 shows a graph of the parameterizing function M(qa). Function M(a) has global 


os ; m, 
minimum VM, and maximum M, at ao M(a)> >. 


Yo 3 
(1) 





Fig. 3. Parameterizing function M (a): 1 — global maximum level M,, 


m 
3 — global minimum level M,, 2 — level of horizontal asymptote es 
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The change in parameter y, in the intervals (0.M,) (My) (SM) (M,%) determines the 


characteristic intervals of the shape coefficient a , under which the conditions for low friction forces are met. 
For low friction forces, frequency function and damping function have the form: 


, 2k(a? -atl) (b(a?-atl)) 

Oa oe) oe ee 
m, +m,o m, +m,o (37) 
b(a’? -at+l 

x(a) =H =a+l) 
m,+m,o 


For high friction forces, at which the aperiodic motion of the system is realized, the frequency function 1s zero, 
and the damping function has two components: 















w =0; 
24] 2g 4l)) 2k(o2-a41 
5 (a) = 22 a+1) | b(o a atl). 38) 
m, +m,o m, +m,o m, +m, 
bio? -atl) I b(2-a+l)) 2k(02-a4)) 
9,(a) =— oN | 
m, +m,o m, +m, m, +m,o 


On the basis of analytical representations of the frequency function and the damping function, characteristic 
variants and features of the extreme properties of the corresponding functions can be determined, taking into account 
the viscous friction forces. 

Discussion and Conclusion. Of interest are the characteristic variants of frequency functions and damping 
functions depending on the conditions of low viscous friction forces. Variants of value y,, that determine the 


characteristic intervals of the shape coefficient a , under which the conditions of smallness of the viscous friction forces 
are met, are considered. 


1. Let y, €(0,M,). The example y, ~ 0.1 is considered. In this case, the conditions for low friction forces for 
any shape coefficient a € (—00..00) are met. Fig. 4 and 5 show the frequency function (a) and damping function 


o(a) for the mechanical elastic-dissipative system with parameters b=1,m, =3, m, =4,k =3. 





-0.5- 6 
Fig. 4. Frequency functions @ (a) : Fig. 5. Damping function o(a@) : 
(1) and (2) are extreme levels (1) and (2) are extreme levels 


Roots of the equation that is equal to zero of the corresponding determinant 
|4p? + Bp +C|=0, (39) 
are p,=0,+jo,, i=1..4, where , 0.91; o, %1.56; @, =—-1.56; wo, -0.91; 6, =-0.14; 6, =-0.44; 
o, ~—0.44; 
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In Fig. 4, the frequency function reaches extreme values equal to the squares of the frequencies w,° ~ 2.46 and 


o, + 0.82. In Fig. 5, the damping function reaches extreme values that are o, ~—-0.44 and o, ~-0.14. The 
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frequency function and the damping function reach their extreme values when the form coefficients are a,* = —0,65 


and a,*=1.15. 


2.Let y, €(V ra . We consider a mechanical system with parameters b=4, y, ~1.33. Fig. 6 and 7 show 
the corresponding frequency function and damping function. The set of coefficients of forms for which the condition of 
low friction forces is satisfied is: (—00,a,)U(@,,0), where a, are the roots of the equation M(a) =y, . For parameters 
b=4, m =3, m, =4, k =3 the roots of the characteristic equation (39) have real parts representing the dissipation 
coefficients, o, —0.57; o, ~—-1.08; 0, »—2.46; o, ~—-0.57 and imaginary parts representing the frequencies 
o, ¥0.73; 0, *0; @, +0; wo, *—0.73. In Fig. 6, the frequency function has a local minimum @,° ~ 0.53 in the 


interval (—co,a,)U(Q,,0). 





oC 

-4 =) 0 2 +t 

ee 
3 
(3) a 7 
VV 
(1) 
(2) O 
-31 G 
-0.1 
Fig. 6. Frequency function: Fig. 7. Damping function: 
(1) is extreme level at point om (3) is extreme level at point OL; ; (1), (2) are extreme levels at 


point on on the two-valued interval 

On the interval (a,,a,), the frequency function is zero. In turn, in Fig. 7, the damping function in the interval (a,,@, ) 
is double-valued and reaches simultaneously two extreme values o, >—1.08 and 0, »—2.46 at point a, *. In the 
domain (—%,a,)U(a,,0), the damping function is single-valued and has one local extremum o, ~—0.57 at point 
Ons 

3. Let y, € (GM , ). We consider a mechanical system with parameters b=6; m, =3; m, =4; k =3.The 
roots of the characteristic equation (39) have real o, >—0.85; o, »—-0.56; 0, »-4.74; 0, ~—0.85 and imaginary 
parts @, 0.36; @, 0; @, 0; wo, ~—0.36. The conditions for low friction forces are met in the interval (a,,a,), 
where a, * 0.63; a, ~ 2.37. In Fig. 8, the corresponding frequency function is positive only on the interval (a,,a,). 
The local extremum of the frequency function is w,’ ~ 0.13. Outside the interval (a,,a,), the frequency function is 


ZCTO. 


In Fig. 9, the damping function is double-valued in the interval (—«,a,) and reaches simultaneously two extreme 
values at point a, *, which are o, »—4.74 and o, =—0.56. In the interval (a,,0) the damping function is also 
double-valued. In the interval (a,,a,), the damping function is single-valued and has one local extremum o, ~ —0.85 


at point a, *. 
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Fig. 9. Damping function: (1), (3) are extreme levels at point OL, 
Fig. 8. Frequency function: 


fivivent ee ea on the two-valued interval, (2) is extreme level at point OL; 
is extreme level at point O, 


4. Let y, €(M,,0 ). For parameters 6=6.52; m,=3; m, =4; k =3, the characteristic equation (39) has 
only real roots o, = —0.50; o, ~—0.88; 0, ~—-0.96; 6, = —5.26 . The corresponding frequency function and damping 


function are shown in Fig. 10 and 11. The interval of fulfillment of the conditions of smallness of the friction forces 
degenerates into an empty set. 





(2) (1), a 
a) ee ee, 





-6 O 
Fig. 10. Frequency function: Fig. 11. Damping function, formed by two branches: (1), (2) are 
case of degeneracy extreme levels of the “upper” branch at points on and OL; ; (3), (4) are 


extreme levels of the “lower” branch at points a, and a, 


In Fig. 10, the frequency function is zero on the whole number axis. In Fig. 11, the corresponding damping 
function is double-valued on the entire numeric axis and has 4 local extrema o, ~—0.50; 0, ~—0.88; 0, ~—0.96; 


o, ¥—5.26 at points a, and a,. 
a ~ 


re PE 
Be 
(h) 


Fig. 12. Pictograms of topological features of frequency function and damping function graphs: (a)-(d) are pictograms associated 
with graphs of frequency functions shown in Fig. 4, 6, 8, 10, respectively; (e)-(h) are pictograms associated with damping functions 
in Fig. 5, 7, 9, 11, respectively 
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The pictograms shown in Fig. 12, compared to the function graphs in Fig. 4-11, reflect a number of 
topological features of the frequency function and damping function graphs. These features include the shape of the 
graph as a single curve, the presence of bifurcation points on the graph of one curve into two, the presence of two 
nonintersecting curves or a “ring”’. 

Thus, the achieved extreme values of the frequency function and the damping function are related to the 
dynamic characteristics of the mechanical oscillatory system taking into account the friction forces. In particular, the 
extreme values of the constructed frequency function and damping function are related to the proper frequencies and 
dissipative coefficients of damped oscillations. In this case, the issue on the existence of extreme values of the 
frequency function and the damping function that do not coincide with the squares of the proper oscillation frequencies 
requires additional consideration. At the same time, it can be assumed that the forms of frequency functions and 
damping functions that determine the modes of free motions of mechanical oscillatory systems with friction are of 
interest for evaluating a wider range of dynamic properties. 

In terms of practical implementation of the possible control of oscillatory modes of mechanical systems based 
on the connectivity coefficient, there are no fundamental obstacles. For example, possible dynamic state control systems 
may include sensors of vibration amplitudes at control points of a vibrating process machine or vehicle. However, the 
construction of such systems requires detailed consideration of a wide range of features related to the technical object. 

In conclusion, the following points can be noted as the conclusions of the presented studies. 

1. For the mechanical vibrating system considering forces of viscous friction, a method of constructing 
frequency function and damping function that are dependent on the connection form coefficient of the free motion 
coordinates, 1s developed. It is shown that the set of extreme values of the frequency function and the damping function 
displays the proper characteristics of an elastic-dissipative mechanical oscillatory system. 

2. It is shown that the frequency function and damping functions for a mechanical oscillatory two-degree-of- 
freedom system with account for viscous friction, can be represented in two variants determined by the conditions for 
the value of the viscous friction forces for a fixed connection form coefficient; for the conditions of low viscous friction 
forces, the values of the frequency function take positive values, and the damping function has one negative component; 
if the conditions of high viscous friction forces are met, the frequency function takes zero values, and the damping 
function has two negative components. 

3. A method is proposed for constructing possible variants of frequency functions and damping functions for 
various values of system parameters based on a parameterizing function that allows determining the regions of values of 
the connection form coefficient in which the condition of smallness of viscous friction forces is met. A criterion for 
classifying frequency functions and damping functions depending on the topological features of their graphs is 
proposed. 

4. The matrix method for constructing the frequency-damping function for a two-degree-of-freedom system 
can be extended to mechanical oscillatory systems considered in different coordinate systems. 

5. As a physical interpretation of the connection coefficient used in the frequency function and the damping 
function, we can consider the linkage in the form of a gear ratio expressed from the ratio of the amplitudes of the partial 
block coordinate oscillations. The ratio under consideration, along with the static state, can be determined for steady- 
state and damped oscillation modes. Thus, a concept is developed in which the starting point for the analysis of a 
mechanical system 1s the linkage. 
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Infinite plate loaded with normal force moving along a complex path 
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Don State Technical University (Rostov-on-Don, Russian Federation) 





Introduction. A technique of solving the problem on an infinite plate lying on an elastic base and periodically loaded 
with a force that moves along an arbitrary closed trajectory and according to an arbitrary law, is considered. 

Materials and Methods. An original method for solving problems on the elasticity theory for plates loaded with a force 
moving arbitrarily along a closed trajectory of arbitrary shape is considered. The problem on an infinite plate lying on 
an elastic foundation is investigated. The plate is loaded with a normal force moving at a variable speed. The load is 
decomposed into a Fourier series on a time interval whose length is equal to the time of its passage along the trajectory. 
The solution to this problem is realized through a superposition of solutions to the problems corresponding to the load 
defined by the summands of the specified Fourier series. The final problem solution is presented in the form of a 
segment of the Fourier series, each summand of which corresponds to the solution to the problem on the action on an 
infinite plate of the load distributed along a closed trajectory of the force motion. The fundamental solution to the 
vibration equation of an infinite plate lying on an elastic foundation is used to construct these solutions. 

Results. A solution to the problem of an infinite plane, along which a concentrated force moves at a variable speed, 1s 
presented. A smooth closed curve consisting of arcs of circles was considered as a trajectory. The behavior of 
displacements and stresses near the moving force is investigated; and the process of the elastic wave energy propagation 
is also studied. For this purpose, a change in the Umov-Poynting vector is considered. 

Discussion and Conclusions. The results obtained can be used in calculations for road design. The study of the 
propagation of the energy of elastic waves from moving vehicles will provide the assessment of the impact of these 
waves on buildings located near the road. Analysis of the behavior of displacements and stresses near the moving force 


will allow assessing the wear of the road surface. 
Keywords: infinite plate, moving load, arbitrary closed trajectory, variable speed, elastic wave energy. 
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Introduction. In many areas of science and technology, problems related to the propagation of elastic waves 
arise. This work objection is to study laws of the propagation of elastic waves that occur under the action of a moving 
load. Problems of this kind were previously investigated by various authors in a number of papers, which considered a 
variety of problem statements with a moving load and proposed various solution methods. 

Often, when solving such problems, a mobile coordinate system is introduced [1—5], or a quasi-static statement 
of the problem is considered [6—12] to exclude time from the number of independent variables. A number of papers use 


the finite element method [11—13]. Interesting results can be obtained using variational [14-16] or direct techniques 
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[17-19]. Methods based on the application of fundamental solutions to the corresponding differential equations were 
used in [20-22] under solving problems on the elasticity theory about a force moving at a constant speed. 
Problem Setting. Consider a differential equation describing vibration of an infinite plate lying on an elastic 


base under the action of vertical force P: 
2 2 a2 P 
NU +e Ce (1) 


Ed’ 


where U is deflection of the plate; D= — 
120-w’) 


; E is Young's modulus; wu is Poisson's ratio; d is thickness of the 


plate; c* = — ; p 1s density of the plate material; k = a ; K is coefficient of rigidity of the elastic base. 


We will consider a solution for which the energy flow is directed from the sources of elastic wave excitation to 
infinity, and assume that force P moves along a closed trajectory y in an arbitrary way. In addition, we assume that 


P = P(s(t)), where s is the arc coordinate counted from some fixed point of the curve y. Obviously, P = P(s(t)) will be 
a periodic function over ¢ with period T , if s(t) is also a periodic function over t with period T . 


Materials and Methods. Consider the fundamental solution to the equation (1). It can be obtained from the 


expression: 


AW+c°OW+kW = 3x x, O(y — y, )o,(t-T), (2) 


where 6,(t—T) = »y O(t-—t—nT). 


The solution to the equation (2) can be obtained by traditional methods, using the limiting absorption principle, 
and presented as a series: 
W(X,Xo2VoVost -—T= >. W, (X,%pqs Vso, ee , (3) 
k=—00 


where w,(x,X),¥,¥y,,) satisfies the equation: 
2 ee 1;,T 
A’w, —c “@,w, + kw, ep Pg) 


It is known that the solution to the equation (1) can be presented as: 


ree: 


U(X.) = | [[ VW Gx, 192 Dot-DP(X,.¥,.Ddx,dy,dT. 


-T/2 R2 
If the moving force is a single concentrated force, which is described by the function: 
P(s(t)) = O(% — x, (s(Q)O(y — ¥,(s@)) 


the solution in this case will look like: 


T 


UX IN== | Wx (orn (t= de. (4) 


Once defined U(x, y,t) , it is possible to calculate displacements and stresses at any point of the plate. 


: Q2kT . 8 , , ; ea , 
If the expression , = — is large enough, it is required to calculate integrals of rapidly oscillating functions 


(3). For this purpose, the quadrature formula was used [23, 24]: 
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where /, are the lengths of elementary segments into which the interval [a; b] is divided; S(x) 1s approximation of f(x) 


by cubic spline M, = S"(x,). 
Research Results. The calculations were performed for an infinite plate loaded with a normal force that 


moved along a closed curve (Fig. 1). 





Fig. 1. Trajectory of the concentrated force 


During the calculations, it was assumed that the plate thickness was d = 0.25 m, parameter c = 221 m/s, 
young's modulus of the plate material E = 232469 H/m’, Poisson's ratio 1 = 0.36, elastic base compliance coefficient 
K=1.864 m™*. The radii that determine the shape of the force trajectory were assumed to be equal to: R; =5m, 
R,=3m, R3=1m. The formula (3) held 120 terms of the Fourier series, and, when calculating the integral (4), the 
integration interval was divided into 120 equal subintervals. 


The law of force motion along the trajectory was described by the function: 


L-si _T 
s(t) ES OL=DD Lng ase, te[0:T]. 


2sin(lZy 2 
; ; T : ; 
We considered the instant of time ft = e when the moving force was at the same point of the trajectory for any 


T. When the parameter 7 is changed, the speed rate of the concentrated force along the trajectory changes. 

To analyze the stress-strain state of the plate, displacements and stresses were calculated in a rectangular 
coordinate system associated with a moving concentrated force. In this case, the axis t of this system was directed 
tangentially to the trajectory t , and n axis coincided in the direction of the external normal to the area bounded by the 
trajectory n (Fig. 1). 

In this coordinate system, the displacement vector and stress tensor can be presented as: 


U =Ut-t+Un-n+Ww-k , 
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ce St-1t+Sn-iia + Stn-(H +H), 


where k is the normal to the plate. 
Fig. 2 shows the variation of the components of the displacement vector Wt,Utt,Unt, the stress tensor 


Stt, Snt,Stnt (second and third graphs) in the axis t, and the change in the same values along the axis n — 


Wn, Utn,Unn and Stn, Snn, Stnn (first and fourth graphs) at z = hf the speed of the force v = 49.3480 m/s, acceleration 


eo) b) 
along the trajectory w,=0 m/s’ , and normal acceleration Ww, = WA = 487.045 m/s’ (the position of the force on the 
1 


trajectory is marked with a red dot in Fig. 1). 
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Fig. 2. Variation of displacements and stresses (v = 49.3480 m/s) 


When studying the energy propagation of elastic waves, the energy flux vector was calculated: 
P=—(S,Uy + S,yUy iE — (Sip, + SyyUty) J 5 
where s, are components of the stress tensor; u, is time derivative of the coordinates of the displacement vector. 


Fig. 3 shows the propagation of elastic wave energy near a concentrated force whose position on the trajectory 
is indicated with a red dot. The length of the vector shown corresponds to the amount of energy passing through a given 


point in space per unit of time, and the direction of the vector indicates the direction of energy transfer. 
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Fig. 3 Energy flux vector (v = 49.3480 m/s) 
The calculations have shown that with an increase in the force motion speed, the behavior of displacements and 
stresses, as well as the nature of the elastic wave energy propagation, changes slightly. Fig. 4 and 5 show the calculation 


results for speed v = 246.7401 m/s, which exceeds the elastic wave velocity in the plate — c = 221 m/s. The effect of the 
concentrated force speed on the distribution of vertical displacements W is shown in Fig. 6, 7. 
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Fig. 4. Variation of displacements and stresses (v =246.7401 m/s) 
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Fig. 5. Energy flux vector (v = 246.7401 m/s) 
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Fig. 8 and 9 show peak value graphs of displacements and stresses depending on the concentrated force motion 


speed. The position of the force on the trajectory at the time under consideration is marked with a dot in Fig. 3 and 5. 
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Fig. 8. Change in maximum displacements depending on the concentrated force motion speed 
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Fig. 9. Change in maximum voltage values depending on the concentrated force speed 


Fig. 10 shows the elastic wave energy propagation near a moving concentrated force. The tangential 
acceleration of the moving force at this moment was equal w, = 1.5503 m/s’. The calculations are performed for the 
instant of time t=7 with the same law of force motion along the trajectory as in the previous case, hence with the 
same law of change in the speed and acceleration of the force motion. 

At this instant in time, the force was at the point of the trajectory shown in Fig. 10, and its speed was 
zero. Fig. 11 shows the variation of components of the displacement vector Wt,Utt,Unt and the stress tensor 
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An increase in the force acceleration also caused changes in displacements and stresses and affected the nature 
of elastic wave energy propagation. Fig. 12 and 13 show the calculation results for the case w, = —155.0314 m/s”. 
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Fig.10 Energy flux vector (w, = 1.5503 m/s’) 
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Fig. 11. Variation of displacements and stresses (w, =1.5503 m/s”) 
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Fig. 12. Variation of displacements and stresses (w, = 155.0314 m/s’) 
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Fig. 13. Energy flux vector (w, = 155.0314 m/s’) 
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Fig. 14 and 15 show peak value graphs of displacements and stresses depending on the tangential acceleration 
of the concentrated force. 
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Fig. 14. Peak value graph of displacements depending on tangential acceleration 
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Fig. 15. Peak value graph of stress values depending on tangential acceleration 
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Discussion and Conclusions. There is a pronounced dependence of the fields of displacements and stresses on 
the speed and acceleration of the force movement at the parameter variation limits considered above. The nature of 
elastic wave energy propagation also depends significantly on the speed and acceleration. 379 
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Sufficiently large values of the speed and acceleration of the force motion were specially considered to test the 
proposed method under such conditions. The results obtained allow us to conclude that the method is quite stable over a 
wide range of variable parameters. 

The use of the proposed method is quite acceptable for solving more complex problems. For that, it 1s required 
that the differential equations describing them provide the analytical construction of the fundamental solution. This 
method is economical and simple since it uses already known problem decisions to build a solution. 
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Introduction. The design technique for a highly efficient technological process of vibration finishing of parts with small 
grooves and holes is presented. The decision is based on the selection of the granulometric characteristics of the processing 
environments. The cross-sectional shape and geometrical dimensions of burrs on typical parts of radio electronic equipment 
(REE) are analyzed. A generalized burr model has been developed. Methodological principles for the selection of particle 
size characteristics of operating environments are determined. 

Materials and Methods. The new classification and coding according to the constructive and technological principles of 
REE components will provide reasonable selection of the equipment, environments and modes when designing their finish- 
ing process. 

Results. A technique has been developed for selecting the granulometric characteristics of operating environments with the 
account of the major technological problems. Based on the design and technological features of the REE components, the 
dependences are proposed for determining the size and shape of the processing granules. The acceptance criteria for evalu- 
ating the results of vibration processing are determined. It 1s noted that one of the major tasks in the vibration processing of 
parts with small grooves and holes is to provide such in-process time at which burrs are removed, and the roughness and 
other surface parameters meet the technical requirements. In this case, the accuracy of the linear dimensions of the pro- 
cessed surfaces should be considered an indicator of quality. Quantitatively, this criterion is assessed on a specific index 
whose calculation considers the largest actual size before vibration processing, the burr height, the smallest allowable size 
after processing, and the tolerance established by the technical requirements. The process efficiency criterion is defined as 
the ratio of the machinability index to the processing time of a batch of parts or the cycle time per part. The proposed crite- 
rion enables to compare treatment processes under validating the solution to technological problems. 

Discussion and Conclusions. The study results enable to confirm that vibration processing in the organic environment 
contributes to the effective removal of burrs and edge smoothing of small-sized parts of electronic equipment with small 


grooves and holes. 
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Introduction. The recently increased demand for precision machining of parts 1s particularly topical in the ma- 





chine and instrument industry, since in these industries, the annual volume of manufactured parts exceeds hundreds of 
thousands of pieces. It should be emphasized that the details of radio-electronic equipment (REE) have a rather complex 


configuration of the outer contour. Most of them are characterized by such non-technological elements as grooves and 
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small holes, deep holes, blind threaded holes, etc. Finishing of such parts in granular media has proved its effectiveness [1, 
2]. Currently, there are no methods for designing finishing and clearing operations for such parts, which limits the wide- 
spread introduction and further improvement of vibration processing. 

Materials and Methods. Based on the design-engineering analysis of REE parts manufactured at the instrument- 
making plant, classification was developed and their coding was performed according to the design-engineering principles. 
The new classifier will provide selecting equipment, environments and operation modes at the design stage of the part fin- 
ishing process with account for their design-engineering features. 

Finishing and clearing processing plays a special role in providing the quality of REE parts, which is difficult to 
perform under modern production conditions due to the complex structural shape of products. 

One of the major methods of forming blanks in mechanical engineering and instrument making is edge cutting 
machining. However, any mechanical processing based on cutting results in defects such as burrs, sharp edges, etc. They 
occur due to the laws of continuum mechanics, so these disadvantages cannot be excluded even when using modern pro- 
cessing centers and optimum processing modes. As a result, functional, aesthetic, and ergonomic problems may occur dur- 
ing the manufacture and operation of parts [3, 4]. This proves urgency of the problem of deburring under finishing and 
clearing operations of REE parts. 

The research work objective is to improve the methodology for designing highly efficient vibration treatment of 
parts with small grooves and holes based on the selection of granulometric characteristics of flexible working environ- 
ments. 

To achieve this goal, it is necessary to solve the following tasks: 

— to analyze cross-section shapes and geometric characteristics of burrs on typical parts of devices, 

— to develop a generalized burr model, 

— to justify methodological approaches when choosing granulometric characteristics of processing environments. 

To study the geometry and parameters of the cross-section profile of burrs, micro-sections were examined using a 
metallographic inverted microscope equipped with the Thixomet Pro system. As the study has shown, the main shape of the 
burr cross section is a triangle. At the same time, the height of the burr is on average 2-3 times greater than the thickness of 


its base, and the linear thickness sizes do not exceed 0.4 mm (Fig. 1). 





Machine building and machine science 


383 


Advanced Engineering Research 2020. V. 20, no. 4. P. 382-389. ISSN 2687-1653 





c) 
Fig. 1. Cross section of burrs in the studied samples. Sample material: a) BrOC4-3; b) LS-59-1; c) AMg6 


In this case, the choice of characteristics of the processing environment has a decisive influence on the surface 
quality of parts and process performance [5—7]. At enterprises, they use working media whose granules can be classified 
according to their geometric shape and size, material, and size of cutting grains, as well as material of a bundle of grains, 
structure and the production method. 

For ease of designing finishing and clearing vibration treatment, a system classification and coding of granules of 
working media 1s proposed. 

Finishing and clearing of REE parts provides for deburring and surface preparation for coatings. In instrument en- 
gineering, various coatings are used: galvanizing, cadmium-plating, nickelizing, chrome-plating, brassing, palladizing, sil- 
vering, gold plating, passivation, etc. To the surface patches of the electrocontact parts of devices, electrolytic coatings (sil- 
ver, gold) are applied. They are characterized by high electrical and thermal conductivity and chemical resistance under the 
increased humidity conditions [8, 9]. Electroplating requires a certain quality of the workpiece. No burrs or sharp edges are 
allowed. 

As the analysis of structures under finishing and clearing REE parts has shown, environments that meet the fol- 
lowing requirements should be used: 

— high density (at least 1.2 g/cm’) with a low weight of the granule; 

— suitable for processing non-rigid workpieces; 

— suitable for processing workpieces with angled surfaces and areas with limited access to the granule of the pro- 
cessing environment; 

— high wear resistance and ability to keep its shape during processing; 

— uniform structure of the granules. 

Taking into account these requirements and the data of the works [10, 11], it should be recognized that the use of stone 
granular media are most appropriate. When crushing the fruit pips and walnut shells, freeform granules with V-shaped edg- 
es are formed. Thanks to this, the granule becomes an analog of a cutting tool. In this case, different areas of the part are 
available for microcutting. 

Research Results. So, environments made of natural materials are optimal for finishing and clearing REE parts 
with small grooves and holes. The widespread use of this approach is hindered by insufficient knowledge of such environ- 
ments and the unavailability of methods for designing finishing and clearing operations using them. 

As arule, under finishing and clearing vibration processing of REE parts, three technological tasks arise. 

First. If it 1s sufficient to remove the burrs and smooth the edges on the workpieces, there is no need to process the 
internal surfaces of the grooves and holes. This is the most typical task that is solved under processing parts. In general, it 
can be solved through applying granules with size R,,,, greater than the size of the largest hole or groove L, 1.e., Ry > L. 
This will eliminate jamming in them. 


The processing environment contact will be considered as the contact of a single granule until the burr is com- 
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Fig. 2. Geometric layout of granule contact with hole edges 


We take into account the wear coefficient of the granule K,,, which largely determines the consumption of the pro- 
cessing environment, as well as the process performance and the quality of the surfaces. The optimal radius of the working 
medium granule is determined by the expression: 

2 2 
Rey = K, — thay 
where h, — thickness of the base of the formed burr, mm; hy — chamfer size that meets the technical requirements, mm; L 
— linear or diametric value of the processed element, mm. 

Compliance with this condition guarantees the burr removal on the outer surfaces and providing the radius of 
smoothing edges required by the drawing. 

Second. This problem occurs when it is necessary to process and prepare for coating the internal surfaces of 
grooves or holes. In this case, the criterion for selecting the particle size should be the size of the smallest hole or groove 
cmin Of the part. 

To avoid jamming of the working environment particles in the holes and grooves, it is recommended to select their 
size from the ratio: 

Da, = 0.6 = 0.7 Lain 
where D,,, — effective diameter of the processing environment granules, mm. 

However, when selecting the size of the working environment particles using this ratio, it 1s critically important to 
analyze its applicability for processing other holes and grooves of the part that are larger than L,,,;,,. If their size 1s equal to 
2D, or 3D,p, it is possible to jam the particles of the working environment, which will interrupt processing the surface of 


the hole or groove. To prevent this phenomenon, the dimensions of the holes and grooves to be processed L should be with- 
L > 3.2D,p (Fig. 3). 


in the range Lyin < L < 1.7D,,, 2.2D.) < L < 2.7D yp, 











Fig. 3. Dimensions of holes that prevent jamming of the working environment bodies: 
Ar — recommended size, An — unacceptable size 


Third. The primary process task is treatment of the interface points of surfaces. Its solution causes the greatest dif- 
ficulties, since it is hard for granules to reach surfaces that are at an angle to each other. One example is thread. If the tech- 
nical requirements of the part specify the radius of the interface of the surfaces, the size of the granule of the working envi- 
ronment must be equal to it or less: 
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As shown in [12—14], due to violation of this condition, three conditional zones can be formed during vibration 
treatment of angled surfaces: 

— a dead zone (no processing occurs in it) 

—- an open zone, 

— an unstable roughness zone (Fig. 4). 





OW L 


Fig. 4. Untreated areas that occur with the wrong selection of granule sizes 


The dimensions of the dead zone and the zone with unstable roughness can be calculated from the formulas: 


1 
p 180°-q i} 
cos————_ 


2 
180°-«a 


CD <R, 





CB = R,, tan 
For processing hard-to-reach surfaces with organic V-shaped particles, we can introduce the concept of the granule 
permeability coefficient, which is determined from the ratio of the mating angle of the surfaces to the angle of the V-shaped 

edge forming the cutting edge: 
Kup > 3: 
Considering the permeability coefficient when selecting the shape of granules will solve a complicated process 

task associated with finishing treatment of surface interfaces. 

To control the finishing and clearing result, it is required to formulate criteria for the validity of parts or their 
batches. For REE parts, these criteria will be the surface roughness parameters and the stability of the quality characteristics 


of the treated surface (Fig. 5). 
<a) 


E a 
< 


Fig. 5. Scheme for determining the processibillity index: 1 — surface roughness after processing; 2 — actual surface roughness; 
3 — burr; EW — waviness deviation; EF — shape deviation; A,,,, — the largest size before processing with account for the burr height; 


d nom 


A min —- Minimum size after processing; TA — tolerance for the processed size 


The main condition for solving this process problem: 
R, —6> [RiminI, 
where R, — the average statistical reading of the data values taken by the indicators, which correspond to the arithmetic 


mean deviation of the reference surface contour; 6 — the scattering field of the values of qualitative indicators, which cor- 
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respond to the values of metal removal rates under processing; [Rimin | — the value of the quality indicators within toler- 
ance. 

Quality indicators within a batch of parts will be random variables distributed according to probabilistic laws. 
Therefore, the criteria for batch validity (1.e., the optimality of the proposed technology) can be set by the inequality pro- 
posed above. At the same time, you need to take into account the defective ratio that does not exceed the probability 
Pi [15]: 

P:(\Ri; — Ril 2 8) < 
where P; — the probability of i-th quality indicator falling outside the boundaries of the scattering field; D; — the calculat- 
ed value of variance R,. 

One of the main tasks under the vibration treatment of parts with small grooves and holes is to provide such pro- 
cess duration at which burrs are removed, and the roughness and other surface parameters meet the technical requirements. 
The quality indicator is the accuracy of the linear dimensions of the processed surfaces. This criterion is quantified by the 
processibillity index (Fig. 5): 


C = Amax—Amin 
O TA 9 


where A ja, — the largest valid size before vibration processing with account for the burr height; A ,,;, — minimum permis- 
sible dimension after treatment; TA — work dimensional tolerance set by technical requirements (Fig. 5). 
The criterion of process efficiency is defined as the ratio of the processibillity index to the duration of the pro- 


cessing of batch of parts, or cycle time reduced to one part: 


The proposed criterion provides comparing treatment processes when justifying the solution to the problems. 

Discussion and Conclusions. The paper analyzes the cross-section shape and geometric dimensions of burrs of 
typical parts of REE devices. A generalized burr model has been created, and methodological principles for selecting the 
characteristics of granules of working environments have been proposed. 

Methods for selecting granulometric characteristics of processing media depending on the main process tasks are 
described. Dependences are obtained for determining the size and shape of processing granules based on the design and 
technological features of REE parts. The validity criteria for evaluating the vibration treatment results are proposed. 

The study results suggest that finishing and clearing vibration treatment in an environment of crushed pips con- 
tributes to the effective removal of liquids, including burrs, as well as the edge smoothing in the radio-electronic equipment 
parts of complex configuration. 
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Introduction. The technological features of the processing crankshafts by the vibration shock method of surface plastic 
deformation (SPD), which is widely used in the technology of manufacturing machine parts, are considered. The 
research objective is to justify the efficiency of the influence of vibration shock hardening treatment on improving the 
quality and performance of crankshafts (CS). 

Materials and Methods. The methodological studies included the validation of the vibration shock processing flowchart 
and the development of an analytic model for assessing the effect of processing on the change in the macrogeometry 
(warpage) of the CS. 

Results. Flowcharts have been developed for volumetric vibration shock finishing and hardening treatment of CS using 
a vibratory machine with a U-shaped work chamber. Its overall dimensions measure alike or exceed the overall 
dimensions of the CS being processed, and ensure the location of the shaft in such a way that its main axis, coinciding 
with the axis of the main bearing journals, is in the zone of the conditional axis of working mass rotation. The surface 
quality parameters were investigated during their processing on the UVG 4X10 vibratory unit according to the proven 
techniques using the dedicated tooling. It has been established that the vibration shock hardening treatment (ViHT) 
enables, due to plastic deformation of microroughness, to obtain a qualitatively new surface microrelief and to reduce 
its initial roughness, to drastically increase the surface microhardness of the CS main and rod journals; at that, it 
changes the stressed state of their surface layer. A calculated dependence is proposed to assess the total warpage of the 
CS strengthened under the ViHT, and its adequacy is confirmed. It is shown that the warpage of the shaft after the 
ViHT is due to the different tension of the rod and main journals of the CS at the level of K, ~ 0.6. 

Discussion and Conclusions. Vibration shock treatment of CS provides an improvement in the geometric and 
physicomechanical parameters of surfaces of the rod and main journals. As a result of processing all surfaces of the 
shaft, warpage does not exceed the permissible values established by the technical requirements. So, we can conclude 
on the efficiency of the considered method of hardening CS with the aim of increasing their operational properties. 


Keywords: crankshaft, surface plastic deformation, vibration shock method, hardening, surface quality, warpage. 
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Introduction. A set of measures to preserve the accuracy of the CS achieved through forming includes 
improving the methods for obtaining blanks, mechanical processing, as well as the introduction of strengthening 
treatment operations into the general process designed to enhance the service properties of these products. One of the 
most common SPD methods that strengthen chamfers and increase the fatigue resistance of CS is running-in and 
caulking [1—7]. However, the CS strengthening through methods of running-in and caulking chamfers is accompanied 
by their warpage, which causes an increase in the main journal runout, the violation of the initial geometric shape that 
requires additional adjustment. In this regard, the development of new processing techniques of finishing and 
strengthening by SPD methods, which increase the reliability and durability of CS under difficult operating conditions, 
is urgent. Of particular interest in this area is vibration shock treatment, which has received wide practical application in 
the technology of manufacturing machine parts [8—15]. In this regard, the research objective is to validate the 
effectiveness of the vibration-shock hardening impact on the quality and performance of CS. 
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Materials and Methods. To achieve this goal, the following tasks are set: 

— to develop flowcharts for vibration shock hardening of all components of CS; 

— to study the impact of vibration shock treatment on the geometric and physicomechanical parameters of the 
surface of rod and main journals; 

— to propose a computational model for assessing changes in the macrogeometry (warpage) of CS under 
vibration shock treatment and validate it. 

The research was carried out on the UVG 4X10 vibratory unit according to the proven techniques. A special 
device simulating a crankshaft was used to determine the quality parameters of the hardened surface. The measurement 
results were processed using mathematical statistics. The warpage of the studied shafts was determined through 
measuring the radial runout of the journals before and after hardening by an indicator reading in 0.01 mm. At the same 
time, the CS was installed with its main journals on prismatic supports and rotated by hand. 

Research Results. Fig. 1 shows flowcharts for volumetric vibration shock finishing and hardening treatment 
of CS using a vibratory machine with a U-shaped work chamber, whose overall dimensions measure alike or exceed the 
overall dimensions of the CS being processed. During processing, the shaft is installed in such a way that its main axis 
coinciding with the axis of the main bearing journals, is in the zone of the conditional axis of working mass rotation. As 
a result of this arrangement of the shaft, the central zone of low activity is excluded, and the cylindrical surfaces of the 
main and rod journals displaced relative to the main axis of the shaft are processed in the zones of medium and 
maximum pressure. CS is installed in the bed along the guiding groves and is held on the supports attached to it. The 
bed with the shaft is immersed in the chamber with the working medium and fixed on racks attached to the frame. 
Supports allow the shaft to scroll around its axis and evenly strengthen under the dynamic impact of the working 
environment (Fig. 1. a, 5) or through imparting an oncoming or passing rotation to it using an independent drive (Fig. | c). 
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a) b) Cc) 
Fig. 1. Flowcharts for processing CS with rotation under the action of the working medium (a, 5) and using an additional drive (c): 
1 — frame; 2 — springs; 3 — working chamber; 4 — bed; 5 — crankshaft; 6 — working medium; 
7 — vibrator, 8 — auxiliary drive 


To determine the preferred modes of processing CS, geometric and physicomechanical parameters of the 
surface of rod and main journals were studied using cylindrical and annular samples made of 45 and 40X steels with 
initial roughness R,=0.16-0.42 um. Processing was carried out by a working medium consisting of steel balls with a 
diameter of 3—6 mm, according to the basic technological scheme shown in Fig. 1. a, b, using a special mandrel (Fig. 2) 
with different amplitude-frequency characteristics and the process time. 





a) b) 
Fig. 2. Experimental mandrel simulating the crankshaft: a) general view; 
b) installation in the working chamber of the vibration unit 
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It is established that under ViHT, a qualitatively new microrelief 1s formed on the surface of the samples, 
whose arithmetic mean deviation of the profile is less than on the initial surface. The change in the arithmetic mean 
deviation of the profile is critically affected by the amplitude of the working chamber vibrations and the processing 
time. So, for the same processing time period t = 20 min, at frequency f = 25 Hz, under increasing the amplitude from 2 
to 3 mm, mean arithmetic deviation of surface profile is decreased by 1.7 times, and a change in the oscillation 
frequency from 20 to 30 Hz at the amplitude A = 3 mm has reduced arithmetic mean deviation of the surface profile by 
1.1 times. 

Taking into account the capabilities of vibratory machines, the established regularities made it possible to 
validate the amplitude-frequency parameters that provide the most effective power impact on the part under the ViHT 
process at the level of A = 3 mm and f = 25 Hz. Processing samples in these modes for 20 min, as shown in Fig. 3, 
reduces the arithmetic mean deviation of the initial surface profile by 2.6 times or 60 %. Increasing the processing time 
to 40 minutes does not cause a significant change in the parameter under consideration. 





0 20 40 60 


Processing time, min 


Fig. 3. Dependence of arithmetic mean deviation of the surface profile 
on duration of vibration shock treatment 


The efficiency of the selected modes is validated by the study results of physicomechanical characteristics of 
the surface layer after processing with steel balls in the mode: A = 3 mm, f = 25 Hz (Fig. 4, 5). Microhardness was 
evaluated on hardness tester PMT3, the stress state was estimated by the value of residual compressive stresses 
determined by the Davidenkov method through cutting ring samples. 





0 10 20 30 40 50 60 
Processing time, min 


Fig. 4. Effect of processing time on microhardness of the surface layer for the ball material: 
1 — steel 40X; 2 — steel 45 
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Fig. 5. Dependence of residual compressive stresses 0, on ViHT duration of samples installed on the 
bushings (1) and on the shaft (2) 


Analysis of the stress state of samples hardened under ViHT (Fig. 5) has shown that the residual compressive 
stresses of samples fixed on bushings and simulating CS rod journals are 35% higher than those of samples fixed on the 
shaft and simulating main journals. This is due to the difference in the intensity of exposure to the processing medium 
in various areas of the working chamber. 

The next stage of research provided for the development of a computational model for estimating the 
macrogeometry (warpage) variation of CS under vibration shock treatment, and its experimental validation. 
V.N. Emelyanov [1] proposed a graphoanalytical method for analyzing the macro deformation of CS after hardening 
treatment by the SPD method (Fig. 6). 


Fig. 6. Scheme of CS warpage after SPD of chamfers of main and rod journals 


The method consists in determining the displacement of the shaft end 6, in section A—A depending on its 
geometric dimensions: r — the distance from the axis of the main journal (MJ) to the axis of the rod journal (RJ); ; — 
the thickness of the crank web; b — the distance from the MJ axis to the end of the CS; ZL — the distance between the 
MJ axes of one crank; C — the distance from the crank web from the RJ to the MJ axis. 

The offset 6, corresponds to the angle of deviation of the shaft axis from the horizontal whose value depends 
on the angular position of the crank web a, which is rigidly connected to the main and rod journals. The interrelated 
parameters 56, and a depend on the average value of the residual compression stresses 0, in the surface layer after the 
SPD and the depth of their occurrence 6,. In addition, it is shown in [1] that the value of warpage of a multicrank CS is 
mainly determined by the degree of hardening of the MJ and RJ chamfers attached to the end web of the end crank. 

Using the considered method, an expression is obtained for calculating the total warpage of CS reinforced 


through ViHT: 

6 (1-v) b 
Oa Oe Kool Oe (0), (1) 
where E — elasticity modulus of the first kind; v — Poisson's ratio; o,, — dynamic yield strength; K,.2 — coefficient 
of adjustment of the value of the RJ residual compressive stresses depending on their diametrical dimensions and 
deformation parameters of the power impact of a part of the working medium; K, = 1.1-1.5 — coefficient of correction 


of the depth of residual compressive stresses; K,, — coefficient that takes into account the difference between the stress 
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state of the MJ and the stress state of the RJ. According to the results of these studies, this coefficient depends on the 
distance from the walls of the working chamber and is K,, ~ 0.6. 

To validate the model (1), experimental studies were carried out on 5 full-scale CS that underwent complete 
mechanical processing. The crankshaft material was as follows: steel 45, hardness after annealing 180 — 228 HB, rod 
journals of 25 mm diameter subjected to the HFC hardening to a depth of 2—4 mm to hardness of 52-65 HRC, while the 
chamfers of 2—3.2 mm radius with roughness R, = 1.6 um remain without heat treatment. Mechanical characteristics of 
the CS material are: v=0.25; E =2-10°MPa; 0, = 360 MPa. The CS dimensions are: r=37.5 mm, h=27 mm, 
L=254 mm, C= 60 mm. 

CS hardening was carried out according to the basic technological scheme (Fig. 1, a, b) on 2x50 UVG vibratory 
unit with a working chamber volume of 50 dm°. The working environment was a mix of hardened polished balls of 3 — 
6 mm diameter from steel ShKh15 of 60-62 HRC hardness. The processing parameters were: vibration amplitude — 
3 mm; vibration frequency — 25 Hz; processing time — 20 min. 

Warpage of the shafts was determined through measuring the radial runout of the main journals before and 
after hardening by an indicator reading in 0.01 mm. The CS was installed by its main journals on prismatic supports and 
rotated by hand. Fig. 7 shows the experimental and calculated values of the CS radial runout after ViHT in sections of 
the shaft end removed from the end crank web of the end crank at the distances: A—A = 123 mm, B—B = 77 mm, C—C = 
17 mm. From these data, it follows that the deviation of the actual values of CS warpage from the calculated values is 
on average no more than 15—20 %, which allows us to recommend the dependence (1) for calculating CS warpage 
under vibration shock treatment. 
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Fig. 7. Dependence of CS radial runout on the location of the controlled section: 
solid line — experimental values, dotted line — calculated values 
Discussion and Conclusions. Vibration shock treatment improves the geometric and physicomechanical 
parameters of the surfaces of rod and main journals of CS. As a result of processing all surfaces of the crankshaft by this 
method, the amount of warpage does not exceed the permissible values established by the technical requirements. This 
allows us to draw a conclusion about the efficiency of the considered method of strengthening CS to increase their 


operational properties. 
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Introduction. The study results on the single interaction under the centrifugal rotary part machining in the abrasive 
discrete medium are presented. Simultaneously with the numerical simulation, experiments were carried out on a 
centrifugal-rotary unit, and the maximum depth of penetration into the surface of the part, the single-track sizes, metal 
removal in one blow of an abrasive granule, were investigated. The removal of metal from workpieces was investigated 
depending on the processing modes, characteristics of the abrasive particle and the processed material. 

Materials and Methods. The dependences for determining the metal removal from workpieces (steels 45, copper Cu- 
OF, and aluminum alloy D16T) are taken into account depending on the grain size (N3) of abrasive particles. The 
process of a single interaction of an abrasive particle and the workpiece surface 1s considered within the framework of 
the dynamic contact problem of the elasticity theory. The authors have carried out finite element modeling of the 
considered structures in CAE ANSYS package. 

Results. The results of theoretical and experimental studies on the metal removal from workpieces depending on the 
grain size of abrasive particles are presented. The technique of their implementation, the tool and equipment used are 
described. The results of theoretical and experimental studies are compared. Their fine precision is established. 
Abrasive tools and processing modes are selected. 

Discussion and Conclusions. The dependences constructed in the work provide determining the rational values of the 
technological parameters of the centrifugal rotary machining (CRM) process. They can be used under designing the 
CRM processes. Therefore, time and financial resources can be saved to achieve the desired surface quality. 


Keywords: metal removal from the workpiece surface, penetration depth, single interaction, metal removal analysis, 
centrifugal rotary machining. 
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Introduction. In mechanical engineering, grinding technology always provides high surface accuracy, which 





is the last step in surface treatment of parts. Materials with high temperature resistance, hardness, and high strength can 
be processed using grinding technology [1]. To increase the efficiency of the centrifugal rotary machining (CRM) 
process, it is required to optimize the model of friction interaction between abrasive particles and the workpiece surface. 

Machining in a free abrasive environment enables to process parts of various shapes, sizes and materials using 
simple and reliable equipment. Temperature in the processing zone is much lower than under grinding, and the 
machining is accompanied by the supply of process fluid. In [2, 3], the thermal texture analysis and workpiece wear 
forecast were performed using the finite element method and the Gaussian process. The studies presented in [4, 5] have 
shown that when a particle of aluminum oxide 1s sliding, the workpiece surface is destroyed, and the metal is removed 
from the part. This paper deals with the process of a single interaction of a workpiece and an abrasive medium under 


CRM. Simultaneously with numerical modeling, experiments were conducted to study the depth of penetration of the 
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abrasive granule into the surface of the part. The effect of the grit size on the metal removal from the processed parts is 
studied. 

The research to determine the impact of technological parameters of CRM processing on quality and 
productivity is not deep. This hinders a widespread introduction of CRM technology into production. To solve this 
problem, it is required to obtain a theoretical model of the CRM process that provides predicting the results of 
processing at the design stage. 

Installation Description. The CRM operation principle (Fig. 1) is as follows: the workpiece (4) and abrasive 
particles are loaded into the working chamber (1); abrasive particles and the workpiece move in a spiral orbit; the 
rotating bottom (rotor) (2) is connected to the engine; the inner surface of the bottom is covered with a wear-resistant 


material to reduce friction. 





Fig. 1. Centrifugal rotary machining flowchart: 1—cylindrical shell; 2 — rotor; 
3 — abrasive media; 4 — workpieces 


Fig. 2 shows a three-dimensional model of the interaction of an abrasive particle and a workpiece. 

In [6], the process of a single interaction of an abrasive granule and a part was studied. The features of this 
process are presented in [7—10]. It is known from [11] that the temperature in the processing zone is low and does not 
change the structure of the surface layer of the workpiece. This paper studies the contact interaction of an abrasive 
particle and the workpiece surface; the abrasive moves at a speed of vg, the angle of contact with the processed 
workpiece surface is @ = 15°—25° (Fig. 2). 

Investigation of single interaction process. Solving the problem of theoretical modeling of a single 
interaction will allow us to study the effect of technological parameters on the CRM. To create a mathematical model of 
metal removal from the surface of a part, it is required to describe impact of the factors on the shape and size of traces 


of interaction between the granule and the surface. 
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Fig. 2. Model of the abrasive particles and the workpiece: 
1 — abrasive particle; 2 — workpiece 


Many papers were devoted to the process of single interaction [12, 14, 15]. Analysis of the works [12-14] 
shows that the interaction of an abrasive particle and the workpiece surface occurs as follows: when a moving particle 
collides with the surface of the workpiece at a certain angle a, the abrasive particle is affected by resistance force P, 
which consists of tangent P, and normal Py. 

The problem is reduced to the movement of a circular cone at a constant speed v > 0. In [5, 6], the maximum 
depth of penetration under a single interaction is determined. 





Arax = 2V,4.R Si OL Pu (1) 
3K,CO, 
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where p,, — material density; k, — coefficient considering the effect of the workpiece surface roughness on the true 
contact area; c — coefficient of bearing capacity of the contact surface; o, — yield strength of the part material; R — 
radius of the abrasive particle; « — angle of contact of the abrasive particle with the workpiece surface; V,, — 
effective velocity of the moving abrasive particle determined from the formula: 


Vip = KpOR,g >» (2) 


where R,, — effective rotor radius;  — rotor rotational rate; x, .— generalized coefficient of effective velocity. 


Based on the work [4], when an abrasive particle collides with a workpiece with the formation of an elliptical 


processing trace, the values of the semi-axes a, b are found from the formula: 
b=./R?-(R-h,,.). 


a= 5 (clgo — fh... +b, (3) 


where f — particle workpiece surface friction coefficient. 
The destruction of the surface layer occurs due to microcutting. Therefore, when calculating the removal of 
metal, it is sufficient to take into account the number of interactions Np, that lead to microcutting [16]: 
O=Np4, 
where Q — metal removal from the workpiece surface. 


Substituting value V,, from the dependencies above, we get: 





S 
O = PP,wtq aR mpu S,,, >4R° (4) 


2 


O'=P Pwtq upu S,,, < 4R° (5) 


Finite element models of abrasive particles and parts were built in the ANSYS program to explore the 


maximum penetration depth under a single interaction based on the parameters of the CRM technology. 
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Fig. 3. Maximum penetration depth under a single interaction 


Fig. 3 shows the distribution of vertical displacement during impact and sliding of an abrasive particle 
fragment on the workpiece surface calculated in ANSYS with account for the plasticity of the part material (bilinear 
model). 

Experimental Research Methods. The study of the process of removing metal from the workpiece surface 
was carried out on cylindrical samples made of various materials: steel 45, aluminum alloy DI6T, copper Cu-OF, 12 
pieces per sample (Fig. 4). The hardness (HB) was measured on the Brinell hardness tester. The hardness (HB) and 
yield strength of the sample materials are shown in Table 1. 

The following were selected as abrasive media: white prisms PT 15x15 conditionally equated to 25 grit 
(Fig. 5 a); porcelain balls 10 mm in diameter conditionally equated to M60 grit (Fig. 5 b); PT 25x25 conditionally 
equated to 12 grit (Fig. 5 d); abrasive white-green cone with 8 grit, d=30 mm, h=30 mm (Fig. 5 d). 

Abrasive granules and blanks are loaded into the working chamber of the CRM. Processing is performed at the 
rotation speed of w=12 rev/s. Every 30 minutes, the treatment is stopped, the samples are removed from the chamber, 
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Fig. 4. Samples for determining metal removal from the workpiece surface: 
a) steel 45; 5) copper Cu-OF, c) aluminum alloy D16T 
Table 1 


Dimensions and mechanical properties of sample materials 


S le di ions, 
Sample materials Hardness, HB Yield stress o,, MPa 


Steel 45 Z 20x20 190-200 P3840 
Aluminum alloy D16T © 30x50 80-82 P2400 
Copper Cu-OF @ 15x20 56-59 P1800 








d) 
Fig. 5. Abrasive particles used in CRM: a) white prisms PT; 5) porcelain balls; 
c) PT 25x25; d) white-green abrasive cone 


Experimental Results. The dependence of the metal removal process on the abrasive grit size (N3) is studied. 
For a comprehensive test of the theoretical model, the comparison based on the results of processing 12 samples of each 
grade was made. The theoretical calculation results using dependencies (4) and (5) were compared to the experimental 
results. Based on the results of experimental and theoretical studies, graphs are constructed. 

Fig. 6—8 show the dependences of metal removal from the surface of the part on the grain size of abrasive 
particles M60; 8; 12; 25; part materials are steel 45; aluminum alloy D16T, copper Cu-OF; the processing mode — the 
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Fig. 6. The dependence of metal removal on N3, the workpiece materials steel 45 
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Fig. 7. The dependence of metal removal on N3, the workpiece material is D16T 
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Fig. 8. The dependence of metal removal on N3, the workpiece material is Cu-OF 
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Metal removal from the surface of the part (g) 





== Theoretical dependence === Experimental values 


Fig. 9. Comparison of theoretical and experimental dependences of metal removal on N3: 1 — copper Cu-OF; 2 — 
aluminum alloy DI6T; 3 — steel 45 

Discussion and Conclusions. Based on the results of theoretical and experimental studies, the dependences 
presented in Fig. 6—9 were found. After analyzing the research results, the following conclusions can be drawn: 

1. The metal removal Q is proportional to the increase in the grit size of abrasive particles. When processing 
at m@=12 rev/s, for steel 45 at N3=M60, Q=0.169g, at N3=25 Q=0.37g; for aluminum alloy D16T at N3=M60 Q=0.267g, 
at N3=25 Q=0.797g; for copper Cu-OF at N3=M60 Q=1.065g, at N3=25 Q=1.789¢. 

2. The effect of the abrasive medium grain size, the mechanical properties of the material and the processing 
mode, reflects correctly the theoretical model of a single interaction. 

3. When comparing the results of experimental studies and theoretical data, the difference is less than 20 % 
(with account for the metal removal on the grain size of abrasive particles N3). 

4. The introduction of research results into production enabled to increase the processing productivity by 1.5—2 
times with the required quality of the treated surface and reduction of the working media wear. 

The results obtained can be used to improve the efficiency of process design of the CRM and to expand 
technological capabilities. 
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Introduction. When administering complex multi-parameter systems, management decisions are often made under un- 
certainty. There is an acute problem of reduction of the likelihood of unwanted events and mitigation of possible dam- 
age. The efficiency of predicting damage to complex systems depends directly on the quality of processing methods, 
systematization, and the amount of input data. It 1s required to improve methods for assessing and predicting damage 
and to develop new approaches and criteria for statistical forecasting of damage and evaluating the system reliability. 
The solution to such problems is complicated by a large number of indicators, data uncertainty, short series of observa- 
tions, incomplete initial information, insufficiently developed scientific methodological apparatus. Existing methods for 
predicting damage in the systems of potentially dangerous objects do not take into account the causes of accidents that 
happened due to unfavorable circumstances. As a consequence, management decisions are made upon unreliable fore- 
casting results. In this regard, an urgent scientific task is the development of methods and techniques for the formation 
of viable management decisions, free from this shortcoming. The major study objective is to consider a particular prob- 
lem for predicting damage due to unfavorable circumstances associated with the indistinguishability of the initial data. 
The tasks are to consider this kind of uncertainty which includes indistinguishability of the true system condition and 
the real value of its quantitative characteristics; to formulate a combinatorial problem for the case when a rather danger- 
ous composite feature is determined by the joint manifestation of two or more simple features. 

Materials and Methods. Under the conditions of multiple indistinguishability, the following was used as the source da- 
ta: a set of indistinguishable outcomes with reliable information on the event instance and the uncertainty of assigning 
the event to a certain type; a family of sets having the same number of elements. The Cartesian product of the families 
of the corresponding sets and the actual value of the group of a compound potentially dangerous factor with a com- 
pound rather dangerous feature are taken into account. The resulting mono-element fuzzy group is presented, which is 
also a possible event resulting from the intersection of two necessary events. 

Results. It is established that the problem of predicting damage due to unfavorable circumstances corresponds to a com- 
binatorial-type problem, which consists in enumerating all sets of arguments. The resulting range, which is an elemental 
group of indistinguishability, characterizes the smaller and larger possible values of the size of the group of a potentially 
dangerous factor with a composite rather dangerous feature. It is shown that the formulated combinatorial problems 
without significant changes are applicable to problems in a generalized form, when composite rather dangerous features 
are determined using not only the operation of intersection, but also uniting and difference; thereby, the initial groups 
are not necessarily the objects with simple features. 

Discussion and Conclusions. The results obtained are focused on the construction of analytical algorithms for establish- 
ing indistinguishability under the monitoring, modeling, forecasting state-related processes and complex dynamic multi- 
parameter objects. 


Keywords: indistinguishability, probability, mathematical model, risk, random event, accident potential. 
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Introduction. Activities in the technogenic sphere involve risks and critical situations. This can be the destruc- 
tion of systems (engineering), loss of control (military matters), bankruptcy (economy). The prediction problem is a 
challenge for many industries and is directly related to the need to improve, develop, and apply the mathematical appa- 
ratus of control tools for complex multiparametric systems. 

Issues of simultaneous processing of dynamic arrays of varying degrees of granularity remain urgent. The pro- 
totypes of mathematical models containing such structures are also optimization problems of practical resource alloca- 
tion under the conditions of possible hard-to-formalize impacts [1]. A high degree of process uncertainty reduces the 
feasibility of using resource-intensive distribution algorithms. At the same time, it is necessary to obtain many alterna- 
tive solutions. This is particularly important in situations when information on existing threats is uncertain and incon- 
sistent [2]. If an object or process changes unexpectedly, an adequate mathematical model may not be available at the 
time of decision-making. At the same time, the failure of complex production facilities often causes man-made emer- 
gencies with severe economic, environmental and social consequences, which makes it necessary to improve the math- 
ematical foundations of risk analysis [3, 4]. 

Strategic decisions on managing complex multiparameter systems are made under conditions of uncertainty. 
The purpose of modern risk management is to avoid a critical situation. In the case of a negative scenario, it is required 
to minimize losses. To predict damage and make adequate decisions, constant monitoring is carried out, which identi- 
fies the major factors of the implementation of the critical situation [5]. 

It is necessary to rank and systematize risks according to the degree of impact on the protected object activity. 
Any system is characterized by many parameters, which are often random and changeable. The efficiency of predicting 
damage to complex systems depends directly on the quality of processing methods, systematization, and the amount of 
input data. Therefore, to estimate the time and condition of the occurrence of a critical situation, combined forecasting 
methods are used, which contain expert, analytical, and simulation parts and use the apparatus of probability theory [6, 
7]. Critical values can be determined analytically based on the experimental data processing results. This is often per- 
formed under the conditions of lack of information, so errors are possible. 

A quantitative characteristic of hazard (or security) of systems and situations is the risk. Hereinafter, the risk 
will be considered damage (consequences) from the implementation of a possible event in a complex system [8]. Due to 
the complexity and heterogeneity of the factors that affect the system and should be taken into account in the decision- 
making process, it is required to use a set of methods for analyzing and processing information to assess the damage 
caused by system failures. 

Often, due to the complexity and high cost of monitoring, some of the parameters associated with the technical 
condition of the system (hereinafter — TCS) are indistinguishable. In general cases, indistinguishability is considered as 
the state uncertainty, which can manifest itself not only in the future, but also in the present (or near past) time [9-10]. 
Indistinguishability is described in terms of the theory of possibilities, in which methods for obtaining simple numerical 
characteristics — estimates of various TCS are developed. It also corresponds to the decision- making principle based 
on facts. According to this principle, the indistinguishability of the system states is caused by the uncertainty of 
knowledge on it and represents the concretization of the knowledge uncertainty on the operation, physical behavior. 
Uncertainties can be exogenous and endogenous, due to external and internal actions, respectively [11]. Each type of 
uncertainty can significantly impair the decision accuracy. 

Current state of the problem. A significant number of papers are devoted to certain aspects of the problem of 
forecasting and risk assessment in case of data indistinguishability. However, these studies lack a general methodologi- 
cal basis. In the papers [12—27], particular solutions are proposed, but the general mathematical apparatus is not devel- 
oped yet. Thus, in [14—15], it is assumed that the evaluation of each of the indistinguishable criteria is equal to the 
arithmetic mean of their numbers. The number of the entire group as a whole object in the ordering is taken for the rank 
of each of the indistinguishable criteria. 

In the paper [16], it is offered to unify approaches to management of complex security of various systems; and 


also, the situation when the expert does not distinguish some criteria is shown. 
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In [17], the problem of constructing interval estimates for an unknown probability in the presence of multiple 
indistinguishable outcomes in the experimental results is considered. Two solution paths are proposed: to take into ac- 
count all unobservable outcomes or to discard them. In both cases, roughening up the result, and errors are possible. 

The paper [18] clarifies the uncertainty and indistinguishability that arise when diagnosing the state of power 
plants. Indistinguishability 1s understood as the uncertainty of the state of the managed object for the observer. At the 
same time, it 1s assumed that uncertainty decreases as we move down the hierarchical ranks of energy system manage- 
ment. In addition, it is proposed to introduce a certain threshold and consider the solutions indistinguishable if the 
square of the difference of the desired value does not exceed the set threshold. 

The indistinguishability aspect is mentioned in light of the development of approximate set theory in [19]. The 
concepts of lower, upper approximation and border area are revealed. This provides creating decision rules “if ... then” 
and using facts only, without assumptions. Referring to [20], the author clarifies that one of the key decision-making 
methods in the field of multi-criterion optimization has been developed on the basis of the theory of approximate sets. 

In the paper [21] devoted to the safe operation of lifting cranes, a method of expert assessment of the frequency 
of adverse events 1s proposed, which provides the development of recommendations for reducing the risk. 

In [22], the problems on military operational research are considered in the following context: warships seek to 
contain and prevent sea robbery, and the model of pirate movement is based on forecasting the probability of piracy and 
on the Markov assumption. The author uses the solution to the flow problem with minimal costs. The number of search 
engines does not matter. It is assumed that they are identical and indistinguishable. 

In the study [23], the T-indistinguishability operators stated geometrically are studied as a special case of gen- 
eralized metric spaces for further application under researching fuzzy subgroups. 

It is shown in [24] that the existing indistinguishability of data from a set of observations does not provide ac- 
curate estimates of the system state. Hence, an accurate prediction should be based on the probability density of indis- 
tinguishable states. This density can be calculated as follows: first, through calculating the maximum likelihood esti- 
mate of the state, and then — through an ensemble estimate of the density of states that are indistinguishable from the 
maximum likelihood state. 

In [25], a characteristic of functions that provide combining partial T-indistinguishability (relations) operators 
into a new set, is presented. In [26], an aggregation of partial T-indistinguishability operators and partial pseudometric 
is considered. The aggregation of a set of partial T-distinguishability operators is analyzed, and the relationship between 
functions is shown. These relations: 

— combine partial T-indistinguishability operators, 

— preserve partial T-pseudometric under the aggregation process. 

In [27], the connection of indistinguishability and fuzzy subsets is validated. It is proved that the basic relation 
for them is lattice isomorphism. 

So, the known approaches to the problem of indistinguishability of influencing factors are functionally limited. 
They do not provide the required level of reliability when making management decisions. The weak points of their re- 
search and methodological procedure are obvious. This generally hinders the development of forecasting systems taking 
into account the factor of indistinguishability, and proves the demand for the development and improvement of the 
mathematical apparatus. 

Thus, an important task of risk management is to predict the damage caused by adverse circumstances due to 
the indistinguishability of the initial data. Its solution is of considerable theoretical and practical interest for many com- 
plex parametric dynamical systems. 

Existing methods of predicting damage in the systems of potentially hazardous objects do not take into account 
the causes of accidents that occurred due to adverse circumstances. As a result, management decisions are made on the 
basis of unreliable forecasting results. In this regard, it is an urgent scientific task to develop methods and techniques for 
the formulation of appropriate management decisions that are devoid of this disadvantage. 

Research Problem Statement. Damage caused by the hazards and attacks is not always possible to predict in 


advance due to the difficulty of distinguishing between the usual states of multiparametric objects and systems whose 
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behavior is not determined. As a result, there 1s a need to study and evaluate two independent types of uncertainty: in- 
distinguishability and nondeterminism. Consider a prediction model due to an unfavorable set of circumstances associ- 
ated with the indistinguishability of the source data. 

We show that in terms of set-theoretic modeling, a mathematical model of adverse circumstances represents a 
combinatorial problem with common source data. We formulate this problem for the case when a very dangerous com- 
pound feature is determined by the joint manifestation of two or more simple features. 

Materials and Methods. Source data: 

— a Set of indistinguishable outcomes U with the number of elements n containing reliable information on the 


fact of the event implementation and the uncertainty of assigning the event to a certain type; 


— family A of sets of the type A; , having the same number n, of elements. 


Under the conditions of multiple indistinguishability, we are forced to consider all sets of type A: , for each i, 
because any of them can represent a group of objects that have acquired the i-th simple feature. 


Consider the Cartesian product of families 4, x A, x...x A, . Its elements are all element sets (4,,4,, ...,4,) of 
sets of type 4, 4,, ...,4,. 

Let us consider the case when the set U = {1,2,3,4} , A=2 of an element set of type A, contains n, =3 ele- 
ments, and of type A, —nlo n, =2 elements. 

Assume, A,, ={1,2,3}, 4,, ={1,4}. Then 

— the true value of the group of a compound potentially hazardous factor with a very dangerous compound 


feature. In this case, the families 4, and A, have the form: 


A, ={{1,2,3},{L2,.4}413,4},{23.4}}. (1) 


A, = {{1,2},{13}, {14},{23},{24}. {3 ay}. (2) 
The composition of groups of a potentially hazardous factor with simple features has the form (3, 4), which 


represents mono multiple indistinguishability groups: 
Ay :{1,2,3,4}, (3) 
Ayy:41,2,3,4}. (4) 
As can be seen, the obtained families (1) and (2) contain C} =4,C; =6 elements, represent multielement 
groups of indistinguishability, and are necessary events. 


The Cartesian product of these families will represent a collection containing the following n,, =C,-C; = 24 


distinct pairs of elements 4, u A, : 


({1,2,3} {1.2} ),({1 2,31 {1,3}),...5({1,2,3} (3,44), 
({1, 2,4} {1,2}),-..5({L2,4} {L.4})>---0( {13,4} {3.4}), (5) 


242) a oxtereaiees: ,({2,3,4} {3,4}). 

You can see that 4, x A, is a family of distinct pairs of sets of 4,, A, , type that have the same number of ele- 
ments — 3 and 2, respectively. At k >2 the Cartesian product A, x A, x...x A, is a family of sets, each of which con- 
tains one set of type 4,, A,,...,4,. 

We will call the mentioned sets — elements A, x A, x...x 4, — k-sets. Each k-set will correspond to the inter- 


section of the sets included in it, that is, we select common elements in such sets. For example, the intersections for the 
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sixth and ninth pairs in (5), which have the form ({1,2,3} {3,4} )and({1,2, 4} {1,4}), will be the sets 
{1,2,3' 0{3,4} = {3 and {12,4 {1,4} = {1,4}, respectively. 


The general result (5) will be: 





{1,2} {1,3} {1} {2,3} {2} {3} 
{1,2} {1} {L4} {2} {2.4} {4} (6) 
{1} {3} {4} {3} {La} 34} | 


12} 13h 14F 123} 124E 1345 
The expression (6) is the resulting mono-element group of indistinguishability, as well as a possible event ob- 


tained as a result of the intersection of two necessary events (1) and (2). 


Now consider function f with arguments that are the mentioned intersections of sets of type 


A, AA, O...0A,, and with values equal to the number of elements of these intersections: 
f(A 94,20...04,)=|4.0.4,0..04,[=9, (7) 


where & — the number of elements of the corresponding set. 


For the sixth and ninth pairs discussed above: 
f({L 2,3} 0{3.4}) = £({3}) <1 
f ({1.204} ft} = f({le}) =2 
Now you can see that we are talking about a combinatorial type problem, which consists in enumerating all 


sets-arguments to obtain different values from them. Such an enumeration can be found after determining the smaller 


and larger values f , that coincide with the smaller and larger possible values of the population of a group of objects 


with a compound potentially dangerous factor. 

Research Results. Consider the following definition: C — a possible event under the implementation of F 
event, ifENC + @. 

We will proceed from the fact that the concept of the possibility of an event is usually associated with assump- 
tions under uncertain circumstances. Therefore, the application of this concept is inappropriate if E, C are known, and 
the realized outcome x € E. With the appearance of £, it is clear that C either takes place, if x € C, or it does not — oth- 
erwise. At the same time, there is no uncertainty. 

The most common, productive, and sufficient concept for our purposes is the concept of a possible event, due 
to the uncertainty in the form of indistinguishability of the outcome that has appeared among other outcomes E. In this 
regard, C is referred to as a possible event under the additional conditions listed below. 

Condition 1. The events E and C are known, and, therefore, their intersection EN C # @ (by definition, this is 
not an empty set). 

Condition 2. The outcome x; is not set, because of which E occurs if more than one element belongs to 
E = {x1,X, ..., X,}. It is only known that x; — one of the elements of £, but it is not established which one exactly, i.e., 
d!x, EE. 

With account for conditions 1 and 2, the generally accepted meaning of the statement that C is an important 
event is specified. So, C either was realized, if indistinguishable actual outcome x; € E/N C (and thus, x; € C), or it was 
not realized — otherwise (if x; € EN C). Consider an example that illustrates a possible event of this kind. 

Let E and C consist, respectively, in the appearance of an even number and more than three points when throw- 
ing a dice, that is E = {2, 4, 6}, C= {4, 5, 6}. Here, EN C = {4,6} + @. It is known that EF occurred, but it was not es- 
tablished which of its outcomes took place. As a result, the above conditions are met, according to which C 1s a possible 
event: it was realized if either the outcome x, = 4, or the outcome x3 = 6, and it was not realized if the outcome x, = 2 


was realized. 
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An important special case of the possible event is necessary event C. It takes place if EN C # @ (condition 1), 
condition 2 is met, and E C C, that is, in the case when all the outcomes of F belong to C. Whichever of them appears, it 
causes the implementation of C. 


In this example, from 7 = 4, n, =3,n, =2, we find that f takes two values: 1 and 2. This indicates that there 
are only one- and two-element intersections 4, ™ A, of sets of the type 4,4, CU, U = tL 2,3,4} , 1.e., In a group of 


potentially dangerous objects with a compound very dangerous feature, there can be either 1 or 2 objects, and the indis- 


tinguishability group (8) is a numerical mono-element indistinguishability group: 
{1,2}. (8) 
The next specific task that follows from the considered model of adverse circumstances, consists in the fact 


that each obtained above value f should be divided by 2, 1.e., finding its normalized values — a set of alternative fre- 


quencies for an arbitrary object to acquire a compound very dangerous feature under the conditions of multiple indistin- 
guishability and nondeterminism. The latter consists in the fact that quotients n;\n represent the frequencies of occur- 
rence on an arbitrary object of the i-th simple feature for all i = 1, ..., k. 


We formulate the purpose of solving the following combinatorial problem corresponding to the frequency 
problem: to determine number s of elements A, x A, x...x A, , to find s(r)/s for all r. 


Quotients s(r)/s for all r represent the probabilities that the number of objects with a compound very dangerous 


feature will be equal to a certain value r. To calculate such probabilities in this regard, we use the dependence: 
GS] CeCe’. (9) 
Similarly, by expanding the scope of a specific problem, it is possible to determine a set of alternative frequen- 
cies of particular interest (absence of a compound very dangerous feature on an arbitrary object). It can be found 


through replacing sets of type 4,, A, with their complements of type U\4,, U\A2 with the number of elements n — 


n, = 1 and n—n, =2, respectively, in the source data. In this case, at the preliminary stage, families, the Cartesian 
product of families whose elements are assigned intersections of sets of the type (U\4,M U\4>), and the function whose 
values are the numbers of elements of such intersections, are defined in the same way. Having determined these num- 
bers, we solve a combinatorial problem corresponding to one specific problem, and after normalization — to another 
one. In terms of modeling the first problem in relation to the new source data, the second one will represent: 

— set U with the number of elements n, 

— various sets U\A4, with the same number of elements n — n,, 

— various sets U\4, with the number of elements n — n> each. 

The combinatorial problem corresponding to the problem with new source data is of great importance for prac- 
tice. To solve it, it is advisable to use the result of the following arguments. Let along with a set of ascending-sorted 
values f, we obtain a set of similarly ordered values — the same as f of function g, whose arguments correspond to the 
new data. It can be shown that the first and second sets have the same number of elements. Let f;, and g, be elements of 
these sets with number h. We can make sure that element g, of the second set corresponds to the same number of argu- 
ments of function g, as element f,, of the first set, and the number of all arguments (variants) fand g are the same, that 
is, equal to s. 

Therefore, the quotients of s(7)/s type, where r = f;, correspond to arbitrary value g, of function g in terms of 
the model. These quotients represent the corresponding probabilities of the presence of groups of objects without a 
compound very dangerous feature of different numbers, if the probabilities of the appearance of different groups of ob- 
jects of the same number with a simple feature 7 for all i = 1...4 are the same. 

It can be shown that the formulated combinatorial problems are applicable without significant changes to prob- 
lems in a generalized form, when compound very dangerous features are determined using not only the intersection 


operation, but also the population 4, UA, U...UA,, differences A,\A,, A, \A,, etc., and the source groups are not 


necessarily only those of objects with simple features. 
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Discussion and Conclusions. Thus, it is established that the problem of predicting damage due to unfavorable 
circumstances corresponds to a combinatorial type problem. It consists in: 

— enumeration of all sets-arguments, 

— derivation of different values /, 

— determination of smaller and larger values f, 

— obtaining the range of possible values of function f, in which its true value is found. 

This range is an element group of indistinguishability, it characterizes the smaller and larger possible value of 
the group size of a potentially hazardous factor with a compound very dangerous feature, and includes several opera- 
tions. 

The results obtained are focused on the construction of analytical algorithms for establishing indistinguishabil- 
ity under the monitoring, modeling and forecasting processes related to the state, and complex dynamic multiparameter 
objects. 
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Introduction. The field of augmented reality (AR) is growing rapidly and has great advances in interaction, navigation 
and tracking. Nowadays there are a lot of trends for AR applications in different areas (education, entertainment, busi- 
ness, medicine, etc.). However, there is a lack of research to provide the evaluating AR apps framework to support de- 
velopers when creating suitable AR applications for specific needs. We provide a practical approach to quantify some 
of the AR applications features. We focus on the development of criteria for evaluating augmented reality applications. 
We discuss the criteria of choosing dimensions for that space such as standards for AR, tools for AR development, nav- 
igation and tracking, content management, usability. We provide analysis and evaluation of AR apps through each char- 
acteristic using guidelines which we have developed. 

Materials and Methods. An AR application is a software application that integrates digital visual, audio and other types 
of content into a real-world environment. The software quality and performance are the main characteristics of the ap- 
plication, which are key factors for AR applications. The analysis of scientific papers, documents and standards made it 
possible to determine characteristics that are the most significant quality indicators based on well-grounded users’ needs 
and demands. 

Results. The criteria we have developed for evaluating applications with augmented reality enable developers to create 
their own software products in stages, based on step-by-step requirements for them, evaluating the development process 
by characteristics. This approach will allow you to create high-quality software products using standardized, modern 
development tools. 

Discussion and Conclusions. In addition, developers will have a detailed understanding of each stage of creating the 
application and the necessary development tools and technologies to obtain the highest quality result. That will give an 
opportunity to decide on specific development tools, methods, models and technologies before starting work on a 
project. As a result, it will provide the final high-quality software product with good extensibility and compliance with 
the modern requirements of the digital industry market. 
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Introduction. The development of the modern market of augmented reality (AR) technology contributes to the 
emergence of a larger number of AR-applications widely used in industries including healthcare, public safety, gas and 
oil, tourism and marketing, entertainment and academia. In connection with the increased interest in this technology, the 
development of functional capabilities of AR projects is also carried out, which stimulates an increase in the needs of 
AR-technology using in the most diverse sectors of modern society. 

Unfortunately, some significant aspects of the AR applications development and implementation and AR ser- 
vices, often ignored for the design simplicity and implementing speed, are the compliance with such systems to real 
conditions and the evaluating under real operating conditions. To develop successful and highly efficient AR systems 
that can be adopted in everyday scenarios, user assessment and feedback are very important [1 ]. 

Augmented reality is a very young industry, and there are still no generally accepted standards for developing 
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AR applications. Although the main platforms for creating augmented reality applications have been defined: ARKit’, 
ARCore’ and Unity® (mostly used) — even now AR applications developed using these tools are available for launch 
only on a limited type of digital and mobile devices. 

The situation is exactly the same with wearable devices, HoloLens and Magic Leap glasses — each manufactur- 
er offers its own unique software for creating augmented reality programs. Experts give several years before the indus- 
try consolidates and common standards allow making the development of augmented reality applications accessible and 
universal for all AR developers. 

Rapid progression of the AR field requires effective and validated methods of design evaluation to be devel- 
oped. Failure to consider the usability of AR applications during the design process will result in an increase in user 
errors and accidents, limiting user trust of the technology and undermining user perceptions of the technology, for both 
AR and Virtual Reality (VR) technologies [2]. 

The authors of the papers [1—3] provided an overview of the important designing and implementing features of 
AR applications and proposed theoretical evaluation of AR systems and frameworks through the standardization as- 
pects. Endsley and others [4] described principles of design heuristics for AR for multi-dimensional augmented envi- 
ronments. Some examples of user experience evaluations were presented in the papers [5—10]. 

J. L. Gabbard and J. E. Swan [11] proposed a Usability Engineering (UE) for Augmented Reality approach that 
inserts iteratively a series of user-based studies into a traditional usability-engineering life cycle. Several usability test- 
ing methods of the AR application (subjective measurement using human perception, objective measure from observa- 
tion, evaluation by expert through cognitive walkthrough, heuristic evaluation, lab observation, questionnaire) were 
described by Pranoto and others [12]. Martins and others [13] presented practical use of the usability methods for evalu- 
ating an AR children’s book with multiple methods. Other aspects of AR technology tools and AR applications evaluat- 
ing were presented in the papers [14] and [15]. 

Akgul and others [16] adapted an existing deep learning architecture to solve the detection problem in AR ap- 
plication using camera-based tracking. Other methods that help to improve AR application and to increment the produc- 
tivity in manufacture were described in [17-19]. 

There are several survey papers on AR development, but none is dedicated to Mobile Augmented Reality. 
Huang and others [20] present the results of the latest technologies and methods survey that improves run-time perfor- 
mance and energy efficiency for the practical implementation of mobile AR applications. 

We have presented AR application “Tilsimli arifler” (“Magic letters”) and special features of design and de- 
veloping mobile AR application for enhancing early literacy skills in our papers [21—22]. 

Software quality and performance. An AR application is a software application that integrates digital visual, 
audio and other types of content into a real-world environment. The software quality and performance are the main 
characteristics of the application, which are key factors for AR applications. Software development is a complex and 
multi-faceted process, in which a large number of specialists of various areas of expertise and various skill levels partic- 
ipate. In addition, many technical, technological, and managerial issues intertwine in the application development pro- 
cess. The success of the project and the quality of the developed product depend on their adequate involvement. 

The analysis of papers, documents (e.g., [2, 3, 10, 16, 22]) and standards (e.g., ISO-IEC JTC 1 SC 24, 18521- 
1*, ISO 9000°) shows that the following characteristics are the key quality indicators based on well-grounded users’ 
needs and demands: 

— inadequate functioning of the software product; 

— insufficient interaction of the product with other software, hardware, telecommunications; 

— failures of the software application during the intended use; 

— the slowed down time of the software product and the delay in the presentation of intermediate and output 
information; 

— incomplete display of information; 

— inconsistency of stored data and information entered by the operator; 

— loss of relevance of the information; 

— violation of the confidentiality of information. 

In addition to such “primary” quality data coming directly from the consumer, developers use “internal” indi- 
cators to evaluate the parameters of the current project: 


' Augmented Reality — Apple Developer. https://developer.apple.com/augmented-reality/ 

* ARCore — Google Developers. https://developers.google.com/ar 

° Unity for all. https://unity.com/ 

* Kim G. Augmented Reality Continuum Concepts and Reference Model — Part 1: ARC Reference Model (Work Item Proposal), ISO-IEC JTC 1 SC 
24, 18521-1. 2012. 

> ISO 9000 Family Quality Management. https://www.iso.org/iso-900 1 -quality-management.html 
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— the lines of code in the standard mode; 

— the number of detected errors per 1000 lines of code; 

— program complexity parameters; 

— the probability of occurrence of specific errors; 

— the project complexity and cost of a code unit; 

— price of a “man-month”; 

— statistical characteristics of processes (expectation, variance, correlation function, etc.) and other estimated 
parameters. 

Augmented Reality Application Evaluation Criteria. Ten experts with more than 3 years’ experience in AR 
application development help us to create criteria for evaluating AR apps. We divided the main criteria for assessing the 
quality of augmented reality applications into 5 groups (Fig 1): 

— AR app design and Art design; 

— Graphic programming; 

— AR applications programming; 

— Application profiling and optimization; 

— Publishing applications (build). 

When developing an AR application design, we recommend the following characteristics that should be con- 
sidered: 

1. The surface. How the application can adapt to various surfaces. If you use frameworks, you can use the 
built-in surface detection function. For example, the application can recognize the floor, walls and objects. 

2. Shine. Evaluation of illumination is very important for the realism of objects. It is advisable to use dynamic 
lighting with shadows in real time. 

3. Space needed. Users can experience AR in three different sizes: table scale, room scale, and open environ- 
ment. The user should always have enough space to simply enjoy the experience. Thus, it is necessary to think over this 
before the user starts using the application. For example, if an application requires an open environment, you need to 
inform the audience in advance before they start using the product. 

4. Single-user or multi-user product. If you are developing a multi-user interface, you should design the prod- 
uct when everyone is involved in the process. It is important to create a sense of audience connection with the product. 
For example, if you are developing an AR game, you can provide a map that shows the location of users and provides 
real-time status notifications. 





Fig. 1. Evaluation criteria for AR applications 


The main criteria are highlighted as follows (Table 1). 
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Table 1 


AR App Design and Art Design evaluation points 


AR App Design and Art Design 


Actual requirements (of the market) in terms of application design 


Interface design requirements 
UX features for AR applications 
Target platform specifications 
3D modeling principles 
Principles of work with particle systems 
Types and features of the creation and use of textures 
Principles of working with shaders and materials, rendering features 
Features for customizing models and textures, as well as materials for export to the game 
engine 
Create high-quality 3D models in modern editors 
Create UV-scan in modern editors 
Texturize models in modern editors 
Rigging and animating models in modern editors 
Customize animation in game engines 
Create and customize shaders, materials in simulation programs and game engines 





The initial version of Table | was proposed by the authors, then it was discussed with ten experts in AR devel- 
opment using brainstorming techniques. 

To solve the problems of displaying complex graphic objects, it is necessary to develop effective methods and 
algorithms for processing graphic information at the stages of input, encoding, transformation and image formation. All 
of this makes up a set of basic computer graphics tasks, so we have identified the main components needed to create a 
high-quality graphic component (Table 2). 

Graphic programming enables to create visual effects through the development of shaders using the functions 
of graphic libraries, customizing the rendering of the development environment, using effective postprocessing librar- 
les, customizing shadows, and more. All these components make it possible to obtain high-quality images by applying 
various effects, as well as improve the optimization of application performance. 

Table 2 
Graphic programming evaluation points 


Graphic programming 


Rendering features on the target platform 
Using textures and materials in the game engine 

Features of the work of graphic libraries 

Principles of rendering geometric objects and images 
Principles of proper postprocessing 
Features of implementing lighting and shadows on the used game engine and target platform 
Configure static and dynamic lighting in the used game engine 
Optimize rendering processes 

Customize postprocessing and final image appearance 
Create procedural geometry using game engine tools 


Programming the rendering of the frame 


Evaluation points presented in Table 2 were proposed by the authors and later they were discussed with ten 





po CReendering features onthe target platform 
| Using textures and materials inthe gameengine | 
po Features ofthe work of graphic libraries 
| Principles of rendering geometric objects andimages | 
p Principles of proper postprocessing 
| __ Features of implementing lighting and shadows on the used game engine and target platform | 3 
| Configure static and dynamic lighting inthe used game engine | 2 
pO ptimize rendering processes 
| Customize postprocessing and final image appearance | 
| Create procedural geometry using game engine tools | 
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experts in AR development. 

During the development of the software part of the mobile application with augmented reality, it is necessary 
to precisely determine the fundamental development tools, and, accordingly, the programming languages, programming 
environments, patterns, architecture, etc. This is necessary for a clearer and more coordinated team work, understanding 


the requirements and tasks, comprehension of the logic and sequence of application development, orientation in the 
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project, reducing risks, and interaction of each team member with each other. To achieve this, we have identified the 
following points for the programming phase, which should be followed in order to obtain a positive result of the devel- 


opment process (Table 3). 


Table 3 
AR Applications Programming evaluation points 


The evaluation points presented in Table 3 were proposed by the authors and brainstormed with ten experts in 





AR development. 

When developing applications with augmented reality, 1t is important to consider the features of application 
optimization for PC and mobile devices, as well as the architecture of mobile devices in the context of application opti- 
mization. In this regard, we have identified the points that we advise to adhere to during the development in order to 
achieve maximum optimization of AR applications (Table 4). 

The profiler tool provides specific data on game performance and facilitates its optimization process. The pro- 


filer provides frame-by-frame metrics with which you can more easily identify problem areas. 


Table 4 
Application Profiling and Optimization evaluation points 


[Oo pimization oF AR application performance id 
Oo piimization oF 3D objets ——SSS~dCi 
Geometry optimization ———SSSS~dSC 


Optimization oF extures and materials forthe target platforms ‘|? 
[Opimization ofthe main application processes | 
[optimization oFphysies in the application ——S—=d 
[Tsing bultsin game engine profilers, as well as extemal profiles [3 





Table 4 was created by the authors and then it was discussed with ten experts in AR development. 

After completion of all stages of the application development, it is important to publish the application. The 
publication represents the release of the application on any platform where the customer (end-user) can easily download 
the final version of product, get acquainted with it, get all the necessary documentations, technical support and feedback 
from the developer. Each of the platforms puts forward its specific requirements to the publishing application, which are 
necessary for correct displaying the application in the platform’s market, obtaining all information about the applica- 
tion’s operation, ensuring end-user security, promoting the application, and more. In this regard, the publication is one 
of the most important and crucial stage of the development. We have identified several main platforms, before using 
which it is necessary to familiarize yourself with all the documentation and assembly features for the appropriate plat- 
form (Table 5). 
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Table 5 
Publishing applications (build) evaluation points 


Publishing applications (build) | 
Features of building an application for Windows 


Features of the build application for Android. 
Features of the build application for OS X/1iOS. 


Evaluation points presented in Table 5 were proposed by the experts in AR development. 





Conclusions. Nowadays, the augmented reality is one of the most innovative and a new digital trend in the de- 
veloping applications for all type of devices. The AR technology opens a new horizon and is going to get more popular 
in the foreseeable future. 

The criteria we have developed for evaluating applications with augmented reality enable developers to create 
their own software products in stages, based on step-by-step requirements for them. This will allow you to create high- 
quality software products using the standardized, modern development tools. 

In addition, developers will have a detailed understanding of each stage of creating the application and the 
necessary development tools and technologies to obtain the highest quality result. That will give an opportunity to de- 
cide on specific development tools, methods, models and technologies before starting work on a project; and, in the 
process of working on the basis of existing criteria, gradually create key application stages, with possible subsequent 
upgrades and improvements. As a result, it will provide the final high-quality software product with good extensibility 
and compliance with the modern requirements of the digital industry market. 

In the follow-up study, we are going to apply this approach for evaluating several AR applications. To thor- 
oughly verify the proposed criteria, additional testing will be required, where more software field experts should be 
involved. Engaging third-party experts will assess the suitability of the proposed criteria. We believe that this approach 


will become part of the AR application development process. 
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Introduction. Verification packages design and analyze the correctness of parallel and distributed systems within the 
framework of various classes of temporal logics of linear and branching time. The paper discusses a polynomially 
realizable class of Ag7 -formulas interpreted on multi-sorted models with hierarchical suspensions. The suspension 
structure is described by an arbitrary context-free (CF) grammar. The predicates and functions of the model signature 
are interpreted on the original CF-list, which is completed during the interpretation process. 

Materials and Methods. A constant model is constructed for theories from Aj) 7-quasiidentities with Noetherian and 
confluence properties. We consider formulas of the multi-sorted first-order predicate calculus (PC) language with 
variables of the “list” sort interpreted on models with a hierarchized suspension. The theory is interpreted in terms of 
grammar inference trees describing the behavior of the specified system. The CF-grammar rules hierarchize the action 
space of the modeled system. It is noted that the expressive capabilities of Ag7-formulas are insufficient for modeling 
real-time systems. Therefore, expressions with unbounded universal quantifier V, known as PT formulas, are used for 
the specification. 

Results. The logical specification of an automated complex which consists of a workpiece manipulator is given as an 
example. The location of the positions is fixed by sensors. The operating cycle of the manipulator is described. The 
specification of its operation consists in the hierarchization of actions by the rules of the CF-grammar and their 
description by the first-order PT-formulas taking into account the time values. 

Discussion and Conclusions. The paper shows that the class of the considered formulas can be used to model real-time 
systems. An example of the logical specification of a manipulator behavior control device is given. 
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Introduction. Mathematically sound, practically significant methods of verification of complex software and 
technical systems are based on the apparatus of mathematical logic [1-4]. The technique of applying this approach to 
various types of real-time reactive systems (communication protocols, control systems, integrated onboard systems of 
space technology, etc.) is known. 

This technique provides verification of model checking systems [5—7]. Numerous verification packages sup- 
port the design and analysis of the correctness of parallel and distributed systems within various classes of linear and 
branching time temporal logics: LTL, CTL, TCTL, etc. [8]. 

To simulate time in these systems, the standard model of the time automaton is used. This is a finite state ma- 
chine equipped with a special type of variable — local clock. Quantitative analysis of the time characteristics of the 
system 1s complicated by complex exponential algorithms for constructing time zones as equivalence classes [9]. There- 


fore, it is required to develop a more expressive, practically significant specification language to simplify the analysis. 
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It is proposed to use the language of X- specifications for simulation, highlighted in the concept of semantic 
programming, which is based on the model-theoretic approach [10]. In this case, the Ist order predicate calculus lan- 
guage, extended by axioms for list structure operations can be used to build a formal model of the analyzed system '*””**. 

Materials and Methods. The paper uses terminology of the papers [11, 12]. Let M be a many-sorted signa- 
ture model o=<I, C , F', R>. Here, J is a set of sorts, including the “list” sort (Jist). C , F , R are sets of constants, 
functions, and predicates, respectively. All signature symbols have a type. If feF' is an n-local function, n => 0, then its 
type 1s <i,, ... i,, i>, where i,,... i, , i¢/, and i,, ... i, are types of arguments, 7 is the type of the function value. Simi- 


larly, n- local predicate reR is of type < i, ... in >. The model carrier M is an indexed family of sets U,=C,, jel, 
where C, is a set of constants of sort 7; f? U, x..xU, 2U,,rc U, x...xU;, . 


For the model M, a list suspension D(C) from the hierarchized CF-lists, whose structure is set by the CF- 
grammar, is formed over a set of constants C. Here, N, T are sets of nonterminal and terminal symbols. The set D,(C) 
is defined as the smallest set of all lists <¢, ,..., t,>, formed for each rule 4 > X,...X, € P, n21 as follows: ¢; is an ar- 


bitrary constant of C, , if.X; €7; otherwise, for X; €N, the element f; is an arbitrary list of X;. 


Ao-formulas are defined in the traditional way as signature formulas ,, using all logical connectives (-, A, V, 
—) and bounded quantifiers Vxet, Axet, Vxct, Ax ct. Here, x 1s a variable of an arbitrary sort; t is a term of 


the /ist sort that does not contain x; ¢€ is the list membership relation; < the nesting relation for lists, defined as 


Ose 3O. > GC Oca. oy IAN 


Below, we will use only bounded quantifiers of the form Vx € y, Vx € y, where nepemennaa y is a variable 


of list sort, the relation € is transitive closure of the relationship €. Denote the indexed sequence of variables x, by x , 


and the membership relation or its transitive closure — by Ee. 

Rules of the CF-grammar hierarchize the action space of the simulated system. For reasons of computational 
efficiency, a class of Aj7-formulas with a “tree” prefix is distinguished. We introduce the relation < — to the “left” for 
the list elements, namely, for the list < --- «,B....> , we consider « ~< f. 


Definition. Aj-formula of the form 


(Vv, ery)... (WVm Elm) (My < by). (tp < I) (V,F) , m= 1, p> 0 
is called Aj7-formula if n;,l; € (v, Tr), 1 < js p; for all prefix variables, the following condition is true: 744 = 7; 
sl sSi<morry, = Uy, k SU. lf 741 = vz, then visp # Vy, and Viz, #7, for all k < i. 

It is easy to show that the prefix of Ay7-formula, due to restrictions on variables, can be presented as a tree 
with root 7,, vertices v,;, 7; and arcs going from vertex 7; to vertex v;,1 <i < m. 

The expressive capabilities of Aj 7-formulas are not sufficient for modeling real-time systems that function cy- 
clically for an indefinite period of time. We will use for the specification of the PT-formula with a universal quantifier 
V. 

Definition. The formula obtained from A)7-formula through ® unbounded universal quantification Vv ®(v), 
is called the PT-formula. 

The model M is defined by a theory of quasi-identities of the form: 


(V7, E 7%)... (V Vm EM) (1, <1, )...(Ny <1, (O71) > W,7)). 
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Biomechanics in Modern University. Rostov-on-Don; Taganrog: SFU; 2015. P. 32. (In Russ.) 

* Glushkova VN. =-specification of real-time robotic systems. In: Proc. Int. Conf. on Algebra and Logic, Theory and Applications. Krasnoyarsk: 
SFU; 2016. P. 22—23. (In Russ.) 

> Glushkova VN. Logical means of diagnosing errors in a hierarchical S-models. In: Proc. Int. Conf. on Algebra and Mathematical Logic. Kazan: 
KFU; 2011. P. 58—59. (In Russ.) 

* Glushkova VN. Logical modeling of robotic technological systems. In: Proc. VIII All-Russian School-Seminar. Rostov-on-Don: SFU; 2013. P. 44. 
(In Russ.) 


Information technology, computer science, and management 


aN 
No 
Uo 


http://vestnik-donstu.ru 


424 


Advanced Engineering Research 2020. V. 20, no. 4. P. 422—429. ISSN 2687-1653 





Here, the formula @ (y) is a conjunction of atomic formulas (or their negations) of the form r, t;=T2, (f=t), f € F,r € 
R, T,, T are terms of signature o. 

The /nt model construction algorithm implements the modus ponens output rule (if @ and gp—>y, then y). The 
input data for the interpreter is a set of initial values of functions and predicates of the form: 
So = {p(C), f(C) = Cnailp € P,f € F}, where C is a set of constants of n elements, n = 1. 

Axioms (ax) are processed in a certain order, first, with positive occurrences of predicates in the left and right 
parts of ax, then, with negative occurrence until fixed points of calculations are obtained. The scope of functions and 
predicates included in the right-hand side of ax, expands when the left-hand side is true. This is because the interpreter 
sets new values for functions and predicates so that the right-hand side of ax is also true. Let the state S,, of the ana- 
lyzed system at the n-th step of the calculation contain the values of all predicates and functions of the model signature, 
and function t P(S) — P(S) in (terminology [3] — predicate converter) reflects the state change when J/nt interpreter 
moves from the n-th step of the calculation to n+1. /nt interpreter constructs the smallest fixed point uZ for the mono- 
tone converter t on P(S). t(Z) =U; t'(Sy), rue T9(Z) = Z, t'*1(Z) = t(t'(Z)) 

Formally, functions f € F and predicate sp € R are interpreted on the CF-list tl(n), which represents the deri- 
vation tree tr(n) in grammar G, where zn is the step of Jnt work. Due to the difficulty of presenting a compact form of 
the interpretation algorithm on the elements of the CF-list, we first give a verbal explanation of the algorithm, focusing 
on tr(7) tree. The input data for the model construction algorithm (/nf) are as follows: the derivation source tree tr(0) 
in grammar G, which is expanded under the construction of model M and Fact=S,). The CF-grammar is used in the 
process of building the model as follows. First, the rules P hierarchize the space of actions and states of the analyzed 
system. We assume that the action names represented in the model signature by predicates and the names of the corre- 
sponding nonterminal grammar symbols are the same. Secondly, the symbols from alphabet V of the grammar uniquely 
define the sorts of all elements of the model universe, including lists, which are assigned to the sort defined by the root 
mark of the corresponding tree. Sorts will be designated mnemonically with initial lowercase characters for the names 
of nonterminal and terminal grammar symbols with the addition of s (sort) symbol at the end. The main advantage of 
CF-grammars is the possibility of using effective syntactically oriented (SO) methods for analyzing the correctness 
(verification) of the model developed in the theory of syntactic analysis of programming languages. 

The interpreter starts by viewing tree tr(0) from the root top to bottom, from left to right. The prefix of all ax- 
ioms satisfies the constraints A)T-formulas. Prefix sorts are defined by the symbols of the CF-grammar G. The inter- 
preter selects as constants the truth domains of the predicates included in axiom ax € Th, the constants associated with 
the vertices of the tree bush, viewing it from top to bottom, from left to right. Moreover, the bush root is marked with 
the name of the corresponding predicate. To reflect the dependence of the simulated technical system on the sequence 
of input signals, it is required to complete the source tree tr(0). To this end, the sequence of rules pr*(ax) € P* is 
attributed to the tree output obtained in the previous step of the algorithm. Moreover, constants from the truth domain of 
predicate r are used as terminal symbols subordinate to the tree vertex marked with nonterminal symbol r. 

We describe the interpretation algorithm Int more formally, without detailing the procedure Con(Q,Th) — 
obtaining all logical consequences from the set of formulas Q based on the axioms of the theory Th. Thy = Sy. Theory 
Thyos S Th includes only positive occurrences of predicates. Thye, S Th includes the negative occurrence of predi- 
cates on the right side of the axioms. We denote t1%,;, t%eg — the derivation trees generated during the interpretation 
process. 

Q:=@; 

Q':=Tho; 

while Q # Q'do 

Qpos *= Q; 

Qnos = Q'; 

while Qyos # Qnos do 
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Qpos *= Qpos: 
Qnos = CON(Qnos, TApos) 
end while 
return (Qpos »tTos) 5 
Qneg = @; 
Qneg = Qpos: 
while Qneg # Qneg do 
Qneg = Qneg; 
Qneg = CON(Qneg: TRneg) 
end while 
return (Qneg: Treg); 
Q = Qpos: 
Q' = Qneg 
end while 
return (Qneg, t™eg) 

The verification of model M consists in checking the properties that the analyzed system should satisfy. We 
express these properties by arbitrary AyT-formulas. Using SO-methods of checking formulas, a proof can be construct- 
ed in the same way as in [13]. 

Theorem. Arbitrary AjT- formula with m-bounded generality quantifiers is tested on the CF-list of power n in 
time O(n™*?),. 

The list power tl is equal to the cardinality of the set{s |s € tl}. 
The estimate is upper, and it can be lowered to linear if you check the formulas using specific SO-methods of language 
processing 

Research Results. We present a logical specification of an automated complex consisting of a manipulator 
maintaining a processing line (¢/) with two positions: loading and unloading of parts (/d and uld , respectively) [14]. 
Sensors record the location of positions. The manipulator functions cyclically starting from the loading position. 

We present a logical specification of an automated complex consisting of a manipulator maintaining a process 
line (tl) (id and uld, respectively) [14]. 

CYCLE 

1. In the initial position /d to load the part, the manipulator lifts the electric drive in 4 seconds. It compresses 
the automated claws and takes the workpiece (2 sec), lowers the electric drive (4 sec) and moves to the right to the ma- 
chine until the position sensor is triggered t/. 

2. To install the workpiece on the machine at position ¢/, the manipulator raises the electric drive, unclenches 
the automated claws (2 sec), lowers the electric drive. Next, the manipulator waits for 4 min, after which it repeats the 
same procedures as in position /d. Then, the manipulator moves to the left to the unloading position until the limit 
switch is triggered uld. 

3. In 8 seconds, the part is unloaded on the conveyor. The manipulator moves to the left until the sensor is 
locked to the loading position /d. Further, the process of complex operation is cyclically repeated. 

The system specification consists of several levels. The manipulator behavior is determined by the signals of 
sensors that record its position: /d, uld, tl (—ld, —uld, —tl, negation indicates the absence of the corresponding sign). 
This sequence of signals is represented by the tuple mc=<x, y, z...>, where x = Id (—/d), y = uld (— uld), z=tl (= tl). It 1s 
generated by a finite-state machine with an initial state x. This automaton is constructed according to a right-linear 
grammar with the rules: 

St Id St, | IdSt)| ¢ 

St; 1 St»|Ad St, 

St.—> uld St| = uld St. 
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Denote a set of lists, made up of symbol strings generated by this grammar, by Dsig. 

The external discrete time (variable n in the logical specification) is determined by the number of transitions in 
the automaton. To describe the second level of the manipulator operation, the CF-grammar is used, which indicates the 
sequence of actions (Oper) of the manipulator: 

— I, La — loading the workpiece by the manipulator in position /d and ¢t/, respectively; 
— Unl, Unla — unloading of the part in position u/d and tl; 

— Mover, Movel — movement of the manipulator to the right and to the left; 

— Lstop, Astop, Ulstop — stop of the manipulator in the corresponding position; 

— Exp — waiting; 

— Cr — failure of the manipulator control device. 

The states of the manipulator (symbol Pos) are affected by its actions. In this example, the state is character- 
ized by continuous time Timec and discrete Timed, given by a natural number. The value of sort Timec 1s segments of 
the form <¢,, t, >, t,, t, — constants, and < is replaced by ( or [ depending on whether the left border is included in 
the time segment or not, similarly for >. 

When specifying the manipulator behavior, you disregard the value of the time of manipulator movement from 
one position to another (determined by signals from position sensors — input to the manipulator control device). The 
signals that are sent to the manipulator actuators for the movement and operation of the automated claws are output sig- 
nals. 

Below are the grammar rules G. 

. Start {Oper>*. 

. Oper» L| La| Unl| Unla| Mover| Movel|Lstop| Astop | Ulstop | Exp| Cr. 
.L- St. 

. La St. 


l 
2 
3 
4 
5. Unl > State. 
6. Unla > State. 
7. Mover — State. 
8. Movel > State. 
9. Lstop — State. 

10. Astop — State. 

11. Ulstop— State. 

12. Exp — State. 

13. Cr > State. 

14. State Timec Timed | Timed Timed. 

15. Timec—> Timed | (Timed, Timed) | | Timed, Timed) | (Timed, Timed | [ Timed, Timed]. 

Timed — a class of tokens whose values are natural numbers calculated under the interpretation of theory 7v. 

In theory Th, variables in formulas are designated mnemonically according to their sort: p (State) = p (State), p 
(oper) = p (Oper), p (n) = p (t) = p (Timed), p (ct) = p (Timec). Predicates Ld, Tl, Unid are defined on the set Timed. 
Ld (n) is true if the manipulator is in the loading position. Similarly, for 7/7 (n) — at the position of the processing ma- 
chine, Uld (n) — at the position of unloading. Let us list the areas of definition of the remaining predicates: Lstop, As- 
top, Ulstop, Mover, Movel, Cr < TimedxTimed; L, La, Unl, Unla, ExpcTimec x Timed. The formulas use the standard 
functions head (<x,,...x, >) = x,, tail(< x,,..x, >) =< x,,...x, >and function Mc: Timed— Dc. Here, Dc is a set of lists 
of sensor signals. 


At the initial time ¢ = 0, n = 1 and at the Ist step of the calculation, predicate Lstop(0,1) is executed; Mc (1) = 


mc, where mceDc. In axioms 1-11, variables t¢, 1 and ct are bound by restricted quantifier V state < oper, V t,n, ct 
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€ state. In axioms 12-17, variable n is bound by restricted quantifier V;_s, = <<<<0,1>>>> — the initial value of the 


list on which the theory is interpreted. Its list constituents, in order of nesting depth, have the following sorts: p (R), 


e (Oper), p (Lstop), p (0) = p (1) = p (Timed). Tree T, corresponds to list s, (Fig. 1). 





Fig. 1. Tree 7) corresponding to list 5, 


In axioms of the theory, the sequence of grammar rules G, that complete tree 7, , is given in square brackets 
on the right. 

Axioms of the theory 

1. Lstop(t, n) \ Ld (n) > L([t, t+ 7), n) A Mc (n + 1) = tail(Mc (n)) [1; 2.1; 3; 14.1; 15.3]. 

2. L(ct,n) A Tl(n + 1) > Mover (ct [2], n + 1) A Astop (ct [2], n + 1) [15 2.5; 7; 14.231; 2.8; 10; 14.2]. 

3. Astop (t,n)— Unla ({t, t + 7], n) [1; 2.4; 6; 14.1; 15.2]. 

4. Unla (ct, n) > Exp ((ct [2], ct [2] + 180), 2) [1; 2.10; 11; 14.1; 15.2] 

5. Exp (ct, n) > La ({ct [2], ct [2] +3), n) A Mc (n + 1) = tail(Mc (n)) [15 2.2; 4; 14.1; 15.2]. 

6. La (ct, n) A Uld (n + 1) > Movel (ct [2], n + 1) A Ulstop (ct [2], n + 1) [1; 2.6; 8; 14.2; 1; 2.9; 11; 14.2]. 

7. Ulstop (t, n ) > Unl ([t, t + 7), n) A Mc (n + 1) = tail(Mc (n)) [1; 2.3; 5; 14.1; 15.2]. 

8. Unl (ct, n) A Ld (n + 1) > Movel (ct [2], n + 1) A Lstop (ct [2], n + 1) [1; 2.6; 8; 14.2; 1; 2.7; 9; 14.2]. 

9. Unl (ct,n) A a Ld (n + 1) > Cr (et [2], n + 1) [15 2.11; 13; 14.2]. 

10. La (ct, n) A A Unld (n + 1) Cr (ct [2], n + 1) [1; 2.11; 13; 14.2]. 

11. LZ (ct,n) \A TI (n4+1) > Cr (ct [2], n + 1) [1; 2.11; 13; 14.2]. 

12. head (Mc (n)) = Id > Ld (n). 

13. head (Mc (n)) = = Id — = Ld (n). 

14. head (Mc (n)) = tl > TI (n). 

15. head (Mc (n)) = =tl > = TI (n). 

16. head (Mc (n)) = uld — Unld (n) 

17. head (Mc (n)) = suld > = Uld (n). 

Theory 7h has the Noetherian property, since the change of the variable under the quantifier V is limited to k 
— the number of elements in the source list mc. In this case, head (Mc (k+1)) 19s undefined, because Mc (k+1)=< >. 
Note that the strings —/d and others, in the right part of axioms 12—17 have the sort string, and — is considered not as a 
logical operation, but as a symbol. 
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For the initial value of function Mc (1) = <l/d, tl, uld, Id, tl, uld> , we obtain a set of consequences: 
Lstop(0,1), Ld (1), L (0,7), 1), Me (2) = < tl, uld, Id, — tl, uld>, Mover (7,2), Astop (7,2), Unla ([{7, 14], 2), Exp({14, 
194], 2), La({194, 197), 2), Mc (3) = <uld, Id, — tl, uld>, Movel (197, 3), Ulstop (197, 3), Unl ({197, 204], 3), Mc (4) = 
</d, — tl, uld >, Movel (204, 4), Lstop (204, 4), L ((204, 211), 4), Mc (5) =< 7. tl, uld>, Cr (211, 5). 

The resulting set of consequences is hierarchized according to the inference in the grammar G obtained as a re- 
sult of the rules assigned to the interpreted axioms. According to them, 12 more vertices, marked with the same symbol 
and connected by edges to the root, are added to tree 7, to the right of the node marked with symbol Oper. Subtrees 


with roots marked with symbols Ld, Mover, Astop, etc., with their states and constants of sort p (Timed) obtained as a 
result of interpretation, are added to the new vertices. 

Discussions and Conclusions. On the constructed model, you can check the truth of arbitrary A, 7-formulas. 
For example, we formalize the statement: “If the manipulator was in the loading position at the instant of time ¢ at step n 
of its operation cycle, then at step n + 2 after197 sec, it starts unloading for 7 sec”. The formula below is tested on a 


given list Oper of sort p (Oper): 
(V state € oper) (VteState) (Vnestate) (Lstop(t, n) > Unl ([t + 197, t + 204], n + 2). 
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Introduction. Electronic document management systems (EDMS) are used to store, process and transmit large amounts 
of information. Automation of these processes is a challenge that requires a comprehensive solution. Its solution will 
reduce the time and material costs for design and make the transition to a more advanced, adaptive EDMS. The paper is 
devoted to the development of new methods for automating the process of distributing information in the EDMS. The 
work objective is to improve the accuracy of the information distribution in the EDMS through moving from analytical 
or algorithmic solutions to the use of new methods based on machine learning technologies. The application of neural 
networks in the furtherance of this purpose will also improve the efficiency of software development through 
automating the analysis and processing of information. 

Materials and Methods. A new method of the automated information distribution based on machine learning 
technologies including a mathematical description of the information distribution rules is proposed. The formulated list 
of conditions for the information distribution provides the implementation of software based on neural networks for 
solving the problem of automatic data distribution in the EDMS. 

Results. The method of automated information distribution has been tested on the example of the EDMS subject area 
when solving the problem of analyzing the correctness of information entered by the user. In the course of experimental 
studies, it was found that the proposed method, based on machine learning technologies, provides better accuracy (8 % 
higher) and is more efficient (in accordance with the Jilb metrics and cyclomatic complexity). 

Discussion and Conclusions. The results obtained confirm the efficiency and accuracy of the method proposed. The 
presented results can be used to automate the processes of distribution and verification of information in adaptive 
EDMS, as well as in other information systems. Based on the method developed, it is also possible to solve connected 
problems: search for duplicates and similar documents, classification and placement by file categories. 


Keywords: electronic document management systems, information distribution, automation, adaptability, machine 
learning. 
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Introduction. Electronic document management systems (EDMS) are widely used for storing, processing and 
transmitting information of various types: documents, spreadsheets, graphic and audio information, engineering and 
accounting documentation, etc. The main functionality for working with documents over the years of the EDMS 
development has been formed. Further development is aimed at automating human activities, increasing flexibility and 
reliability, and implementing intelligent decision support modules [1-3]. 

However, such flexibility and adaptability of the EDMS causes additional time and material costs, increases 
the complexity of designing and upgrading the system. Therefore, the automation of system design is a challenge for the 
development of adaptive information systems, including EDMS. This implies many separate tasks for automating the 
processes of analysis, processing and transmission of documents, the solution of which in totality provides reducing the 
load on developers during the implementation of adaptive functions in EDMS. 


Obukhov A. D. Automation of information distribution in adaptive electronic document management systems using machine learning 





Within the framework of study, the issue on classification of information and its subsequent automatic 
distribution into specified categories in the EDMS will be considered. This sub-task is one of the most common when 
organizing electronic repositories and archives, filing cabinets, and filling out document forms. Distribution refers to the 
placement of data and files in specified positions: by categories, media, directories, and so on [4]. Classical algorithms 
are not always able to detect an error in the placement of information without involvement of an expert or moderator. 

Considering approaches to solving the problem of information distribution, we can distinguish approaches 
based on the use of machine learning’. The application of artificial intelligence technologies for solving classification 
problems has proved its efficiency in numerous studies and experiments [5-8]. When solving the problem and 
implementing the information distribution method, machine learning technologies will be used to automate the data 
classification process. 

Based on the analysis of the information movement process, the following urgent tasks can be identified for the 
classification and distribution of data in adaptive EDMS: 

— classification and determination of compliance of the entered data with the category in which the user 
placed them [9]; 

— data categorization’ , and, in case of an error, moving them to the correct categories, or raising a warning 

to the user [10]; 

— definition of data duplication by features [11]. 

To solve these problems, it is proposed to develop a method for automated information distribution in adaptive 
EDMS, which generalizes existing approaches to information classification and is based on the use of machine learning 
methods for automating data processing. 

Materials and Methods. We formalize the main stages of the method of automated data distribution in 
adaptive EDMS. 

The method is based on the formation and training of a neural network for information classification. 
Therefore, at the first stage of the method implementation, it is required to prepare a set of information objects 
X = {x,,...X,} for traming and testing the neural network. The object x, ¢ X can be represented by text or numerical 
information entered by the user in the form fields, or by files uploaded via the EDMS interface [12]. The data 
preparation process can include normalization, tokenization, lemmatization, and extraction of file properties and 
attributes. For the collected data, a set of categories Y = {y,,...,y,,} 18 predetermined, that is, there is a continuous 
display X — Y . We approximate this mapping by neural network NN : 

NN(X)=Y. 

With a sufficient amount of training data, it is possible to provide required accuracy of the classification. In 
accordance with the studies presented in [13], at least 50 copies of training data should be provided for each output 
feature. Then for Vx, e X , we get NN(x,)=y,, where j=1.M. 

At the second stage of the method, information 1s distributed. Let the EDMS have a set of information objects 
X = {X,,....X,} and a set of corresponding categories Y = {y,,...,y,,} . To distribute information, the following 
conditions should be checked: 

— if Vx,eX, NN(x,)=y, and object x, are placed by the user in category y,, then the position x, remains 
unchanged; 

— if dx, «eX, NN(x,)=y,, (where y,, — the category of malicious objects), then object x, is removed from 
the EDMS, a warning is sent to the system administrator, the user who added object x, , is blocked; 

— if dx, eX, NN(x,)=y,, but object x, is entered by the user in the category y, 

(y, #y; ), then it is required to redistribute the information in the EDMS. 


Consider all possible options for redistribution: 


1. Category y, is free (empty), then object x, 1s transferred from category y, to y,. 


' Umadevi S, Marseline KSJ. A survey on data mining classification algorithms. In: 2017 International Conference on Signal Processing and 
Communication (ICSPC): IEEE, 2017. P. 264-268. 

* Popova ES, Spitsyn VG, Ivanova YuA. Using artificial neural networks to solve the problem of text classification. In: Proc. 29th Int. Conf. on 
Computer Graphics and Vision “Graphicon’”’. Bryansk, 2019. P. 270-273. (In Russ.) 
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2. Category y, is occupied by some object x, (such that NN(x,)=y,), then object x, is transferred to the 
buffer, the user is given a warning about incorrect input x, and duplication of information. 

3. Category y, is occupied by some object x, (such that NN(x,)#y,), then object x, is transferred to the 
buffer, and object x, is transferred from the category y, to y,, the user is given a warning about incorrect input x, . 


The stages and conditions of information distribution considered in the framework of the method include the 
main scenarios for adding information to the EDMS. Through the use of machine learning technologies, it 1s possible to 
automate the distribution of information into specified categories. 

Research Results. To implement and integrate the method into adaptive EDMS, we will use a microservice 
approach. The neural network is implemented using the Keras library (Python), then it 1s imported into a microservice, 
which can be implemented on the basis of any framework, for example, Flask. Data between the EDMS and the 
microservice is transmitted through the HTTP Protocol in JSON format, which provides compatibility with any EDMS 
implementation [14]. 

As an example of the data distribution task, we will use a document card form in the EDMS with 5 fields: 
document name; document author; contact information (address) of the author; date of creation; document description. 
In the course of the pilot study, the input of user data into the form and filling in the information in the specified 5 fields 
will be simulated. The neural network will classify the entered information and check it for compliance with the 
specified category. 

To create a set of source data, we use the generator based on open database of Russian names, cities and 
countries, as well as sequences of Russian words and figures. This will enable to get constructions that correspond to 
real data, but do not store personal data of real people. Text data will be processed using a tokenizer to convert it to a 
numeric format. The maximum sentence length is limited to 20 words. We generate a training sample. We will generate 
10,000 elements for each category. We will also add 50,000 items with incorrect and erroneous data. Thus, we get an 
array of 100 thousand elements. 

The relationships between input attributes and categories are shown on the heatmap in Fig. 1. The map 
provides visualization of the correlation matrix and performing a visual analysis of the data to determine how the 
parameters affect each other and the output variables [15]. 
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Fig. 1. Raw data heatmap 
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The network learning process is shown in Fig. 2. The final accuracy of the neural network on the test set after 
5 epochs amounted to 97.8%. 


1.00 


Accuracy 


mms ACCuracy, 


1.0 1.5 2.0 25 3.0 3.5 4.0 4.5 5.0 
Epoch 
Fig. 2. Classification neural network learning process 


The number of errors did not exceed 7 % for the “Author” category and 1 % for the “Description” category. In 
other categories, incorrect data was recognized in 100% of cases. 
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Fig. 3. Testing the neural network for incorrect data entry 
The second experiment (Fig. 4) consists in entering data from ; categories to i ( j #i). A set of 100 elements 


for each category i consisting of elements of other categories (25 for each category) is formed. The graph shows the 
number of tests in which the neural network incorrectly recognized data from other categories as corresponding to the 
current category. In the second experiment, 1 % of errors was allowed in the “Name” category. The remaining 
categories were worked out without errors. 

Next, we compare the accuracy of the proposed method (hereinafter referred to as the “neural network 
method’’) with the classical solution to the problem of classification and distribution of information (referred to as the 
“classical method”). For the classical method, the following results were obtained: 16 % of errors in the “Name” 
category, and 25 % — in the “Description” which is much worse than indicators of the neural network method (Fig. 5). 


On average, the use of neural network methods provides accuracy increase by 8 %. 
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Fig. 5. Comparison of accuracy of neural network and classical methods 


Next, we will compare the methods by metrics of software implementation complexity. The cyclomatic 
complexity expressed by McCabe number (total and averaged) [16], as well as the complexity according to the Jilb 
metrics [17], are used as metrics. The final results are presented in Table 1. 


Table 1 
Comparison of classical and neural network methods 


Classical method Neural network method 


Cyclomatic complexity (averaged) A (3.71) A (1.37) 
siitbmewies ——]SSC~—~“ RSC“‘S(SS!C*'YOOC~C‘“‘CS#CONSSCOOCS™S 





Thus, the complexity of the neural network method in terms of cyclomatic complexity and Jilb metrics is 
lower. The accuracy of the developed method is almost 8% higher. It is worth noting that the complexity of 
implementing the classical and neural network methods of data classification and distribution is relatively small, so, 
they can be considered comparable. However, the obtained results on classification accuracy confirm the high 
efficiency of the neural network approach. 
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Discussion and Conclusions. The paper sets the task of automating the processes of data classification and 
distribution in adaptive EDMS. The proposed method of automated data distribution includes the formalization of the 
classification and distribution of data, the use of machine learning technologies to automate the solution to the problem. 
The formulated list of information distribution conditions provides implementing the software based on neural networks 
that solves the problem of automatic data distribution. 

To test the developed method, experimental studies were conducted on the basis of generated data on 
documents in the EDMS. The accuracy of the trained neural network was about 98 %. Additional tests have shown that 
the neural network can detect incorrectly distributed data in almost 100 % of cases, and, under the worst conditions, the 
error did not exceed 7 %. Thus, the efficiency and accuracy of the proposed method is confirmed. The developed 
method, in comparison to the classical implementation based on algorithmic support, shows the following positive 
effect: an increase in average accuracy by 8 %, a decrease in the complexity of the implementation. 

The presented results can be used to automate the processes of distribution and verification of information in 
adaptive EDMS, as well as in other information systems. In addition, on the basis of the proposed method, it is possible 
to solve related tasks: search for duplicates and similar documents, classification and placement by file categories. 
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Introduction. This work is aimed at solving the problem of phytoplankton dynamics in the coastal environments using 
the example of the Azov Sea. This takes into account the transformation of forms of phosphorus, nitrogen and silicon, 
as well as the aquatic medium motion, the distribution of temperatures and salinities over the sea area. River flow, 
varying in volume and chemical composition, affects significantly the variability of hydrophysical and biogeochemical 
parameters of the processes occurring in the coastal environment. This explains the need for statistical processing of the 
data from long-term observations over the river flow characteristics. 

Materials and Methods. The mathematical model of biogeochemical cycles is based on a system of non-stationary 
equations of the convection—diffusion—reaction of parabolic type with nonlinear functions of sources and lower-order 
derivatives, to which the corresponding initial and boundary conditions are added. In the course of statistical analysis of 
the series of long-term observations over river flows, the values of the following indicators were found: skewness 
coefficient, degree of kurtosis, variance and standard deviation, coefficient of variation, autocorrelation coefficient, 
Neumann ratio, and Anderson criterion. 

Results. The statistical analysis of the series of long-term observations over the hydrochemical indicators of the Don 
river suggests heterogeneity of the field data. This is due to the stochasticity of nutrient inputs and the volume of 
freshwater flow to the sea as a result of natural and anthropogenic factors. Field data should be correlated with seasonal 
changes in the aquatic environment temperature. This paper presents the results of a computational experiment to model 
the dynamics of phytoplankton populations in summer season, when temperatures are favorable for their reproduction 
and growth. The proposed mathematical model considers the spatially inhomogeneous distribution and transformation 
of forms of phosphorus, nitrogen, and silicon, as well as changes in salinity, temperature, and motion of the aquatic 
environment. 

Discussion and Conclusions. The multispecies mathematical model of the dynamics of phytoplankton populations is 
considered with account for the transformation of forms of phosphorus, nitrogen, and silicon in the coastal 
environments. The analysis of data from field observations, for which its major statistical parameters are calculated, is 
carried out. As a result, it is concluded that data of the long-term observations are significantly variable. This is due to 
two reasons. Random nature of the input of nutrients and the volume of river flow as a result of anthropogenic factors is 
the first reason. The second reason includes the alternation of relatively high-water and low-water periods for fresh flow 
over the last 12-15 years. The hydrological regime is changing mainly due to the reduction of the average annual 
freshwater flow of the Don and partly of the Kuban. This trend is likely to increase due to climate changes, as well as 
with further regulation of the Don river flow after the Bagaevsky hydroelectric installation start-up. Numerical 
experiments based on the field data confirmed the predictive validity of the developed models and programs. They can 
be used to predict change in the composition and abundance (concentrations) in the Azov sea core planktonic 
populations, which define, on the one hand, food resources, and, on the other hand, the aquatic environment in terms of 


the ongoing sea salinization. 


Keywords: biogeochemical cycles, phytoplankton population, biogenic substance, chemical-biological source, 


convection—diffusion—reaction equation, field data. 
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Introduction. The sea of Azov is a large estuarine coastal system. It is the shallowest sea in the world. It 
warms up almost evenly in summer (with temperature differences on average no more than 4 °C). At the same time, it is 
characterized by a large difference in salinity — from 0% to 12-15%, since river flows provide an influx of fresh water 
commensurate with the total volume of sea water, and salty Black sea water comes from the Black sea near the Kerch 
Strait. River runoff affects significantly the biochemical composition of the water body 

[1]. Mathematical modeling of biogeochemical processes, which provides performing diagnostic and 
prognostic calculations of the marine ecosystem dynamics, is considered important. A variable volume and 
hydrochemical composition of the river runoff affects significantly the parameters of hydrophysical and biological 
processes occurring in the coastal system. Therefore, it is advisable to conduct a statistical analysis of the data of long- 
term observations, in particular hydrochemical indicators of runoff of the rivers running into the sea of Azov, and to 
predict biogeochemical processes on the basis of statistically processed input data. 

In the area of hydrodynamics and forecasting of marine systems, the works of Marchuk G. I. [2], 
Matishov G. G., Sukhinov A. I. [3], Berdnikov S. V., Tyutyunov Yu. V. [4], Yakushev E. V. [5], Ilyichev V. G., and 
others should be mentioned. The paper presents the results of complexing a mathematical model of biogeochemical 
cycles and a model of hydrodynamics of the sea of Azov [6-8]. This provides improving the modeling accuracy and 
considering such factors as hydrodynamic processes in coastal systems, heterogeneous distribution of temperatures, 
salinities and biogenic substances that affect the development of phytoplankton populations, the transition of biogens 
from one form to another [9]. It should be noted that numerical models of spatial-three-dimensional hydrophysical 
processes in coastal systems are the subject of a separate study of the authoring team. They provide considering the 
dynamically changing geometry of the bottom and coastline, wind stress on the free surface and its elevation, friction on 
the bottom, Coriolis force, turbulent exchange, evaporation, river runoff, deviation of pressure values in the aquatic 
environment from the hydrostatic approximation, etc. In this paper, the input data (distribution of the three-dimensional 
velocity vector, as well as salinity and temperature) are the results of numerical calculations based on the hydrophysical 
model [10]. 

Materials and Methods. To describe the model, an initial-boundary value problem is formulated for a system 
of parabolic equations with lower derivatives and nonlinear right-hand functions: 

OFi 4 9 SHE yy OH 


Odi _ 7 
r ate v et reiicas div(kgradq;) + R,,, (1) 


where q; = concentration of the i-th component [mg/l]; ic¢M, M= {F, F2, F3, POs, POP, DOP, NO3, NO2, NH4, 
Si}; {u,v,wt — components of the water flow velocity vector [m/s]; 4 — turbulent exchange coefficient [m/s]; Ra. 


function — source of nutrients [mg/(I-s)]. 
In the equation (1), index 1 indicates the type of substance (Table 1). 
Table 1 
Biogenic substances in the model of phytoplankton dynamics 
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Chemical-biological sources are described by the following equations ( i € {1,2,3} , where 1 — is ChV , 
2— AF-A,3— Sc,and ChV, AF — A,Sc — symbolic designations of plankton species): 
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Here, K,, — specific respiratory rate of phytoplankton; K,, — the specific rate of decay of phytoplankton; 
K,, — Specific extinction of phytoplankton; K,, — specific rate of autolysis POP; K,, — phosphatofication 


coefficient POP; K,,, — phosphatofication coefficient DOP; K,, — specific rate of ammonia oxidation to nitrites 


during nitrification; K,, — specific rate of oxidation of nitrites to nitrates in the process of nitrification; s,, s,, S., — 


normalization coefficients between the content of NV, P, Si in organic matter [11—12]. 
The growth rate of phytoplankton is determined from the expressions 


Cu. = Koy Sr (T) fe (S) min} fo (dro,)> Fv (von » Fon» Donn E> 
C,, ~ Kvn Sr (T) fs (S) min \ fp ae 1, ty (Ayo, , D no, Wn, Pe (4s; )} ’ 


where Kyr- — maximum specific growth rate of phytoplankton. 
Dependences of temperature and salinity: 


—_ ‘i ie (s) = > es 
f, (r) = exp| —a r., » f(s) = exp} —B 5 : 


opt 








where 7,,,, S,,, — temperature and salinity, optimal for this type of phytoplankton; 0 > 0, B>0 — coefficient of the 


opt ? O, 


interval width of tolerance of phytoplankton to temperature and salinity, respectively. 





Fig. 1. Model scheme of biogeochemical transformation of phosphorus, nitrogen and silicon forms 


Below are the functions that describe the content of biogens. 
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For phosphorus /, Ce ) = a , where K,,, — constant of half-saturation with phosphates. 
POA POA 
For silicon fs (qs, ) = — , where K, — constant of half-saturation with silicon. 
+Ko 
Si Si 
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where K NO, nitrate half-saturation constant;; K NHy ammonium half-saturation constant; K asi — inhibition 


coefficient of ammonia. 
Assume that the coefficients included in the expressions for the source functions are positive and independent 


of time f. 
For the system (1), an initial-boundary value problem is set in a cylindrical domain G_. Let the boundary 


» of the cylindrical region G be a piecewise smooth surface and Y= X,, UX, Uo , where X,,— the surface of the 
reservoir bottom, &,— the undisturbed surface of the aquatic medium, o — lateral (cylindrical) surface. 
Let u, — be the normal component of the velocity vector of the water flow with respect to 2 , vector of the 


external normal to 2% . For example, for concentrations g, at the lateral border: 





g,=90, ono, if u, <0, ie; (2) 
Og. 
ti _¢), ono, if u, 20,iEM; (3) 
On 
Oq. ; 
—=0, on 2, — water surface, ic M; (4) 
Iz 
aq, | 04; 
ti E1559 VE {P,P 5: eels = €5 iis 
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ic {PO,4, POP, DOP, NO3, NO, NH4, Si} at the bottom ~,,. (5) 


Here, &,,;,&,,; — nonnegative constants; ¢,,, i¢{F'1, F'2, 3} take into account the sinking of algae to the 


ol 


bottom and their flooding; E,,,1€{PO4, POP, DOP, NO3, NO2, NH4, Si} take into account the absorption of nutrients 


by bottom sediments. 
For a system of equations, it is required to specify at any time the velocity vector of the water flow, the salinity 
and temperature field, as well as the initial values of the functions 4g, : 


G(9V20) = Gy 0952). (x,y,z)EG, t=0,ieM, (6) 


V(x, y, 2,0) = V, (x, y,Z), T (x, y,z,0) =T1, (x,y,z) ; S(x, y,z,0) = 58, (x,y,z) ; 

Statistical processing of data from long-term observations on river flows to the sea of Azov. Considerable 
river runoff relative to the volume of the sea affects significantly the biological and hydrophysical processes occurring 
in the sea of Azov [13]. A large amount of biogenic substances, including nitrogen, phosphorus and silicon, the main 
nutrients for phytoplankton, enter the water body from the river flows. In the twentieth century, the main part of the 
water inflow to the sea of Azov falls on the Don runoff — 63 % (Fig. 2-4)’. 
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Fig. 2. Series of long-term observations of the Don river runoff (1993-2012): nitrogen concentration (NV-NH,) 


' Ecological Atlas. Black and Azov seas. OC “Rosneft”; Arctic scientific center; NIR. Moscow, 2019. 464 p. (In Russ.) 
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Fig. 3. Series of long-term observations of the Don river runoff (1993—2012): phosphorus concentration (PO4) 
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Fig. 4. Series of long-term observations of the Don river runoff (1993-2012): 
silicon concentration (SiO 4) 


Table 2 
Results of calculating the statistical parameters of field data 
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The study on the series of long-term observations of the Don runoff allows us to draw a number of 
conclusions. 

— The nutrient concentrations considered have both positive and negative asymmetries. 

— Random variables for nitrogen and silicon are shifted relative to the distribution center, as evidenced by the 
high value of the asymmetry coefficient. 

— Large values of variances and standard deviations were obtained for all biogens. 

— The autocorrelation coefficients are small; therefore, a strong nonlinear trend is characteristic of the full- 
scale data series. 

— The variation in all rows is greater than 20 %, so the rows are highly variable. 

— For nitrogen, the presence of an autocorrelation relationship according to the Anderson criterion is obvious 
for the number of values of 20 in the sample, since the autocorrelation coefficient exceeds 0.299 at a confidence level of 
5%: 

— For nitrogen, the presence of autocorrelation of residues according to the Neumann criterion 1s obvious at a 
confidence level of 5% for 20 observations, since the Neumann ratio is less than 1.2. 
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— Phosphorus and silicon do not show autocorrelation relationships; we reject the hypothesis of 
autocorrelation of residues. 

As a result of statistical analysis of the field data [14], we can conclude that they are highly variable. This is 
due to the stochastic flow of nutrients from the Don runoff and the significantly changing volume of runoff under the 
impact of natural and anthropogenic factors. To use field data in the model (1) - (6), it 1s appropriate to take into 
account seasonal behavior. Hereafter, when modeling, we will consider the summer period. 

Results of numerical experiments. A numerical simulation of the solution to the problem on the dynamics of 
phytoplankton populations in summer with account for the transformation of phosphorus, nitrogen and silicon forms, is 
carried out on the example of the sea of Azov. The simulated area corresponds to the physical dimensions of the sea of 
Azov: length — 355 km, width — 233 km, space step in horizontal directions — 1000 m. Fig. 5 shows a satellite image 
of the sea of Azov, confirming the compliance of the study results with full-scale data. The image shows clearly the 
distribution of green and blue-green algae in the Taganrog Bay area, and diatoms — in the central part of the sea. 





Fig. 5. Satellite image of the sea of Azov taken with a moderate resolution spectroradiometer 
(MODIS) on the NASA Aqua satellite on July 31, 2004 


B pe3yIbTaTe BbIUHCIMTeIbHOrO 93KCIICPHMeHTa TOUyIeHbI CeTOUHbIe paciipeyeieHuA KOHIUeCHTpalui 
OCHOBHBIX MOMyIAWHH PUTOMIAHKTOHA UM MMTATeCJIbHbIX BeIeCTB B A30BCKOM Mope (puc. 6). Ilepuoy pacueta — 30 
CYTOK. DTOFO JOCTATOYHO JIA YCTAHOBJICHHA CTALMOHApHbIX PexXHMOB B 3aa4yax JMHAMMKU (UTOMIAHKTOHA. 
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Fig. 6. Distribution of concentrations: a) green alga Chlorella vulgaris; b) blue-green alga Aphanizomenon flos-aquae; c) diatom 
442 Sceletonema costatum; d) phosphates; e) nitrates; f) dissolved inorganic silicon 


Sukhinovy A. L, et. al. Modeling biogeochemical processes in the Azov Sea using statistically processed data on river flow 





The figures reflect the dynamics of phytoplankton populations, the cycles of phosphorus, nitrogen and silicon. 
In the process of extinction and decay, phytoplankton releases phosphorus in dissolved and weighed organic forms, then 
in the process of phosphatofication, they go inorganic — phosphates, which are consumed by phytoplankton. The 
nitrogen cycle is also described: in the process of vital activity, phytoplankton releases nitrogen in organic form, which 
decomposes to ammonia. Under the nitrification, ammonium is oxidized to nitrites, and then to nitrates. It 1s worth 
noting that phytoplankton consumes all three forms of nitrogen. The consumption and release of silicon by diatoms was 
stated. Comparison to the simulation results for high-water periods shows that in recent low-water years, the habitats of 
green and blue-green algae in the Taganrog Bay area have significantly moved (by many kilometers) eastwards, closer 
to the Don river, a source of fresh water. 

Discussion and Conclusions. The paper presents a multi-species mathematical model of the dynamics of 
phytoplankton populations, form transformations of biogenic substances-phosphorus, nitrogen, and silicon compounds. 
The model takes into account: 

— impact of salinity and temperature on the development of three main types of phytoplankton (green, blue- 
green and diatomic algae); 

— phytoplankton uptake of phosphates and nitrogen forms; 

— transition of forms of phosphorus and nitrogen from one to another; 

— absorption of silicon by diatoms; 

— advective and microturbulent movement of the aquatic environment; 

— water flows and sources at the border. 

To study the field data, a statistical analysis procedure of long-term series of observations of biogenic 
substances concentrations (phosphorus, nitrogen, their compounds, etc.) that enter the sea with the Don runoff has been 
developed and adapted. The statistical analysis of data from long-term observations, in particular hydrochemical 
indicators of the Don runoff, became the basis for predicting biogeochemical processes with account for the movement 
of the aquatic environment, the distribution of temperatures and salinity. The of numerical experiment results are 
consistent with the data of space sensing of the sea of Azov, which confirms the predictive value of the models used and 
the methods of their numerical implementation. The comparison of the distributions of green and blue-green algae 
populations in the Taganrog Bay for high-water and low-water periods shows that in recent low-water years, their 


habitats have significantly moved (by many kilometers) eastwards, closer to the Don river, a source of fresh water. 
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